System and method for pain monitoring using a multidimensional analysis of physiological signals

ABSTRACT

The present invention is for a method and system for pain classification and monitoring optionally in a subject that is an awake, semi-awake or sedated.

CROSS-REFERENCE TO RELATED APPLICATIONS

This application is a Continuation of U.S. application Ser. No.12/779,963, filed May 14, 2010, which is a Continuation-in-Part of U.S.Provisional Application Ser. No. 61/180,161, filed on May 21, 2009 andPCT International Application No. PCTIIL2008/001493, filed on Nov. 13,2008, which claims the benefit of U.S. Provisional Application Ser. No.60/987,782, filed on Nov. 14, 2007, all of which are hereby incorporatedby reference in their entirety.

FIELD OF THE INVENTION

The present invention relates to a system and a method for painmonitoring by performing a multidimensional analysis of a plurality ofphysiological signals, and in particular, to such a system and method inwhich pain monitoring, classification and identification is establishedfor individuals exhibiting various states of consciousness.

BACKGROUND OF THE INVENTION

Pain is an unpleasant sensory and emotional experience associated withactual or potential tissue damage, or described in terms of such damage.The inability to communicate verbally does not negate the possibilitythat an individual is experiencing pain and is in need of appropriatepain-relieving treatment (www.iasp-pain.org/AM/). Pain is alwayssubjective where each individual learns the application of the wordthrough experiences related to injury in early life. Biologistsrecognize that those stimuli which cause pain are liable to damagetissue. Accordingly, pain is that experience we associate with actual orpotential tissue damage. It is unquestionably a sensation in a part orparts of the body, but it is also always unpleasant and therefore alsoan emotional experience. Experiences which resemble pain but are notunpleasant, e.g., pricking, should not be called pain.

“Pain Threshold” is defined as the least experience of pain which asubject can recognize as pain. Traditionally, this threshold has beendefined as the least stimulus intensity at which a subject perceivespain. Properly defined, however, the threshold should be related to theexperience of the patient, whereas the measured intensity of thestimulus is an external event. Because the threshold stimulus can berecognized as such and measured objectively, it has been common usagefor most pain research workers to define the threshold in terms of thestimulus, even though it is preferable to avoid such a definition. Inpsychophysics, a threshold is defined as the level at which 50% ofstimuli are recognized. Thus, the pain threshold would be the level atwhich 50% of stimuli would be recognized as painful. As the stimulus isonly one aspect of pain, it cannot be a measure or a definition of pain.

“Pain Tolerance Level” is defined as the greatest level of pain which asubject is prepared to tolerate. As with pain threshold, the paintolerance level is the subjective experience of the individual. Thestimuli which are normally measured in relation to its production arethe pain tolerance level stimuli and not the level itself. Thus, thesame argument applies to pain tolerance level as to pain threshold, andit should not be defined in terms of the external stimulation as such.

Pain may be described as either a symptom or an indication of anunderlying problem. However pain in and of itself may be a considered adiagnosis or condition. Fibromyalgia is an example of a conditionwherein pain is not a symptom but rather a finding. Many such diagnosisare becoming more prevalent as pain gains recognition for being acondition and not merely a symptom that may be subsided once theunderlying problem is treated.

However, to date state of the art pain monitoring has primarilymanifested and centered on individual in the unconscious state, inproviding a Depth of Anesthesia (herein referred to as DOA) readingand/or monitoring.

Depth of Anesthesia monitoring (herein after referred to as DOA or DOAMonitoring) uses physiological signals that represent certain autonomicnervous system activity or brain activity for monitoring a certain stateof a patient under anesthesia. DOA monitoring is a general term for painand/or awareness and/or muscle activity monitoring when a patient isunder general anesthesia. In the unconscious anesthetized state pain andawareness are difficult to be distinguished as they both may result inthe same physiological symptoms.

Conversely, pain monitoring attempts to detect sensation of physicaldiscomfort and is not limited to the state of consciousness of asubject. Therefore during pain monitoring a patient can be and is oftenfully awake.

Although DOA monitoring has gained in popularity over the last decadeprimarily because of the increase in the number of publication relatingto “awareness during anesthesia”, it is only in the last few years, thatpain monitoring has become a subject to increased awareness.

State of the Art DOA and/or pain monitors are described in U.S. Pat. No.6,117,075 to Barnea, U.S. Pat. No. 6,571,124 to Storm, US PatentPublication No. 2006/0217615 to Huiku, U.S. Pat. No. 6,757,558 to Lange,etc. each describing the independent use of a physiological signals suchas skin conductance, EEG, ECG, PPG, temperature etc., to determine theDOA or pain level. However, medical studies have shown that a usage ofcombination of parameters from different physiological signalssignificantly improved the pain and no-pain classification performanceachieved compared with discrimination using any single signal alone(Guignard 2006),

Other state of the art DOA monitors, such as those described by U.S.Pat. Nos. 6,685,649 and 7,367,949 and European Patent No. EP1495715 toKorhonen describe a DOA monitoring system for a user that is undersedation or anesthesia. These publications are centered and rely uponanalysis of a single parameter associated with the cardiovascularsystem, specifically using blood pressure (BP), heart rate variabilityare correlated to the detection of pain.

U.S. Pat. No. 7,215,994 to Huiku discloses a method for monitoring astate of anesthesia or sedation by comparing cortex related EEGbiopotential signal data from the patient to subcortex-related biosignaldata from the patient, the subcortex-related biosignal data including atleast bioimpedance signal data. However although a few signals are usedtogether to obtain a DOA reading this system is limited to individualsthat are fully sedated and therefore unconscious.

A state of the art pain monitoring system is described in U.S. Pat. No.7,407,485 to Huiki, that presents a pain monitoring system that is basedon one or more physiological parameter that are measured, normalized andthen compared to ‘a threshold surface’, while the frequency of thresholdcrossing infers the relative pain level experienced.

Other prior art publications such as US 2006/0217614 to Takala et al, US2006/0217615 to Huiki et al, US 2006/0217628 Huiki, and in US2007/0010723 Uutela et al, report the use of a group of physiologicalfeatures to form an Index of Nociception to determine the state of apatient.

Prior art teaches DOA systems that are associated with a fewphysiological signals and parameters to infer the pain state of apatient while under sedation. Similarly prior art teaches painmonitoring systems that are limited in that the system is heavilydependent on the number of physiological variables used to obtainappropriate pain indication.

SUMMARY OF THE INVENTION

The present invention overcomes these deficiencies of the background byproviding a system and method for pain monitoring, for pain experiencedin various states of consciousness that optionally and preferably takesinto account a plurality of physiological signals and parameters todetermine the pain level.

Within the context of the present application the following terms areused as is understood and known by those skilled in the art.

States of consciousness within the context of the present applicationthe terms states of consciousness optionally includes sedated, partiallysedate and awake.

Within the context of the present invention the following terms andcorresponding shorthand are interchangeably used throughout the text forthe following terms as is understood and known by those skilled in theart. Standard deviation (‘std’), Very Low Frequency (“VLF”), LowFrequency (“LF”), High Frequency (“HF”).

Within the context of the present invention the term Galvanic SkinResponse (herein referred to as GSR) may also be referred to as ElectroDermal Response (EDR) or Skin Conductance Response (SCR), commonlyrefers to methods for measuring the electrical resistance of the skinand measured, optionally this is measured with two or three surfaceelectrodes and acquiring a base measure.

Within the context of the present invention the term Electro-Gastro-Gram(herein referred to as EGG) is a noninvasive method for the measurementof gastric myoelectrical activity using abdominal surface electrodes.

Within the context of the present invention the term Pupil DiameterMeasurement (herein referred to as PD) measures pupil size and movement.Optionally PD may be measured by infrared videography or computerizedpupillometry.

Within the context of the present invention the term Electromyography(herein referred to as EMG), refers to a technique for recording andevaluating physiologic properties of muscle activity either at rest orwhile contracting. EMG signals are optionally and preferably recordedwith surface electrodes. Optionally and preferably a plurality oflocation specific EMG signals may be recorded from various locations ona subject and/or muscle groups. For example Frontalis (scalp)Electromyogram (herein referred to as FEMG) measures over the frontalismuscle underlying the forehead.

Within the context of the present invention the term PhotoPlethysmoGraph(herein referred to as PPG) is a non-invasive transducer to measure therelative changes of blood volume from a finger.

Within the context of the present invention the term Electro-Cardio-Gram(herein referred to as ECG) is a graphic representation of anelectrocardiograph, which records the electrical activity of the heart.

Within the context of the present invention the termElectroEncephaloGraph (herein referred to as EEG) is the measurement ofelectrical activity produced by the brain as recorded from electrodesplaced on the scalp.

Within the context of the present invention the term ElectroOculaGraph(herein referred to as EOG) is the measurement of electrical activityproduced by the eye movement and retina as recorded from electrodesplaced on the face and frontal lobe.

Within the context of the present invention the term Blood pressure(herein referred to as BP), refers to a signal capture an arterial bloodpressure, i.e., to the force exerted by circulating blood on the wallsof the larger arteries. Optionally BP may be measured by invasive ornon-invasive methods.

Within the context of the present invention the term Laser DopplerVelocimetry (herein referred to as LDV) quantifies blood flow in tissuesfor example such as skin to extract different features such as thevasomotor reflex (SVMR).

Within the context of the present invention the term Capnograph refersto a device provided to monitor the concentration or partial pressure ofcarbon dioxide (CO2) in the respiratory gases in this context it is alsoused for any measurement of concentration end-tidal Nitrous oxide (N2O),oxygen (O2), or Anesthetic Agent.

Within the context of the present invention the term Accelerometer is adevice for measuring acceleration and gravity induced reaction forces.

Within the context of the present invention the term physiologicalsignals refers to any measurable signal or event that is optionallymeasured directly or indirectly from a subject through sensors,transducers or the like preferably providing a measurement indicative ofthe state of a patient. Optionally and preferably the physiologicalsignals may be further analyzed, processed, or otherwise manipulated toprovide further details regarding the state of a patient. Physiologicalsignals for example include but are not limited to, blood pressure,respiration, internal and/or surface temperature, pupil diameter, GSR,and signals received and/or abstracted and/or derived from ECG, PPG,EOG, EGG, EEG, EMG, EGG, LDV, capnograph and accelerometer or anyportion or combination thereof. Preferably a physiological signal mayfurther comprise any signal that is measurable and/or detectable from asubject.

Within the context of the present invention the term feature extractioncommonly refers to the processes, manipulations and signal processingmeasures performed to analyze a physiological signal most preferably toabstract from the signal valuable information, data or subsignalreflective of the state of a patient,

Within the context of the present invention the term feature referrersto at least one or more of physiological features, a priori features, apriori data, external data or external input. Optionally features may bequantitative or qualitative or the like.

Within the context of the present invention the term a priori data, apriori feature, external feature or data may be interchangeably used torefer to any data received or otherwise obtained about a subject forexample including but not limited to nociception response, conceptualresponse, context relevance response, Behavioral response, subjecthistory, gender, type of medicine, diagnostics, patient condition,patient definition of pain level, age group, weight, height, historicaldata, drug history, drug interaction or the like data. Most preferably,this inclusion of a priori data increases the classification efficiencyof the optional classifiers used with the system and method of thepresent invention.

Within the context of this application the term physiological featurespreferably refers to features extracted from a physiological signalthrough feature extraction to ascertain data associated with aphysiological signal. An exemplary list of features optionally utilizedwithin the system and method of the present invention are depicted inTable 1, preferably a plurality of features may be used in anycombination thereof.

TABLE 1 list of optional plurality of features Num- ber of fea- # SignalFeature Description tures  1. PPG PPG Peak (P) and Trough (T) amplitude,mean amplitude and std of amplitude The amplitude of the Peak (P) andthe Trough (T) of the PPG signal, mean amplitudes and STD of amplitudesin predefined time window. Peak denotes a point of maximum blood volumein a finger; Trough denotes a minimum basal blood volume. 6  2. PPG PPGMaximum Rate (MR) point The amplitude of the maximum rate point (MR) ofPPG signal, mean amplitude and STD of amplitude in predefined timewindow. Maximum rate is a point between onset injection and Peak wheremaximum rate of blood volume increase is observed. 3  3. PPG PPGdicrotic notch The amplitude of the dicrotic notch of PPG signal, meanamplitude and STD of amplitude in predefined time window. 3  4. PPGPP/PT/ PN/NT/ NM intervals, mean and std (variability) of interval Thetime interval between peak to peak, peak to trough, peak to notch, notchto trough, notch to maximum rate, in PPG signal, mean interval and STD(variability) of intervals in predefined time window. 15  5. PPG Freq.P-P Variability (PPG-HRV) VLF, LF, MF and HF Power of the Very LowFrequency (0.0033 Hz-0.04 Hz), Low Frequency (0.04 Hz-0.15 Hz), and HighFrequency (0.15 Hz- 0.4 Hz) frequency bands of power spectrum of the P-Pinterval (PPG based Pulse Rate Variability power spectrum) in predefinedtime window. 3  6. PPG Freq. P-P Variability LF/HF Ratio between LF(0.04 Hz- 0.15 Hz) PPG-HRV power spectrum and HF (0.15 Hz-0.4 Hz)PPG-HRV power spectrum in predefined time window. 1  7. PPG Area UnderCurve (AUC) The integral of single beat of PPG signal. 1  8. PPGSpectrum PPG envelope Power Spectrum of the envelope of PPG signal.Envelope—Peak- Trough of PPG signal. 1  9. PPG PPG Variability waveletanalysis Wavelet analysis of the P-P interval variability. 1 10 PPG FreqPPG-RSA (Respiratory sinus arrhythmia) The frequency of dominant peak atHF (0.15 Hz-0.4 Hz) band of PPG-HRV power spectrum in predefined timewindow. 1 11 GSR GSR amplitude, mean amplitude and std of amplitude Theamplitude of the GSR signal, mean amplitude and STD of amplitude inpredefined time window. 2 12 GSR GSR Peak (P) amplitude, mean amplitudeand std of amplitude The amplitude of the Peak (P) of GSR signal, meanamplitude and STD of amplitude in predefined time window. 1 13 GSR PPinterval, mean and std (variability) of interval The time intervalbetween peak to peak of GSR signal, mean interval and STD (variability)of intervals in predefined time window. 1 14 GSR Phasic EDA, amplitude,mean amplitude and std of amplitude The amplitude of the firstderivative of the GSR signal (EDA phasic), mean amplitude and STD ofamplitude in predefined time window. 1 15 GSR Spectrum of the GSR signalPower Spectrum of the GSR signal in predefined time window. 1 16 GSRPeak Amplitude The amplitude of the highest peak of the power spectrumin predefined time window. 1 17 GSR Peak Frequency The frequency of thehighest peak of the power spectrum in predefined time window. 1 18 GSRGSR wavelet analysis Wavelet analysis of the interval. 1 19 ECGQ/R/S/T/P amplitude, mean and std of amplitude The amplitude of theQ/R/S/T/P pulse, mean amplitude and STD of amplitude in predefined timewindow. 15 20 ECG RR/PQ/ PR/QT/ RS/ST interval, mean and std(variability) of interval The interval between each pulse or betweeninternal pulse waves RR/PQ/PR/QT/RS, mean interval and STD (variability)of intervals in predefined time window. 15 21 ECG Freq. R-R Variability(ECG-HRV) VLF, LF, MF and HF Power of the Very Low Frequency (0.0033Hz-0.04 Hz), Low Frequency (0.04 Hz-0.15 Hz), and High Frequency (0.15Hz-0.4 Hz) frequency bands of power spectrum of the R-R interval (ECGbased Heart Rate Variability power spectrum) in predefined time window.4 22 ECG Freq. R-R Variability LF/HF Ratio between LF (0.04 Hz- 0.15 Hz)ECG-HRV power spectrum and HF (0.15 Hz-0.4 Hz) ECG-HRV power spectrum inpredefined time window. 1 23 ECG Freq ECG-RSA (Respiratory sinusarrhythmia) The frequency of dominant peak at HF (0.15 Hz-0.4 Hz) bandof ECG-HRV power spectrum in predefined time window. 1 24 ECG Freq. RRIVariability wavelet analysis Wavelet analysis of the R-R intervalvariability. 1 25 ECG Freq. Alpha Slope of HRV power spectrum inpredefined time window 1 26 ECG Freq. Beta Slope of the log of HRV powerspectrum in predefined time window 1 27 ECG- PPG ECG- PPG PTT PulseTransition time The time interval between R peak of ECG signal and Peakof PPG signal (PTT or rPTT), mean interval and STD (variability) ofintervals in predefined time window. 3 28 Temp Temperature amplitude,mean amplitude and std of amplitude The amplitude of the temperaturesignal, mean amplitude and STD of amplitude in predefined time window. 229 Temp Temp Peak (P) amplitude, mean amplitude and std of amplitude Theamplitude of the Peak (P) of temperature signal, mean amplitude and STDof amplitude in predefined time window. 1 30 Temp PP interval, mean andstd (variability) of interval The time interval between peak to peak oftemperature signal, mean interval and STD (variability) of intervals inpredefined time window. 1 31 Temper- ature Spectrum of the temperaturesignal Power Spectrum of the temperature signal in predefined timewindow. 1 32 Temper- ature Peak Amplitude The amplitude of the highestpeak of the power spectrum in predefined time window. 1 33 Temper- aturePeak Frequency The frequency of the highest peak of the power spectrumin predefined time window. 1 34 Respi- ratory Upper peak amplitude, meanamplitude and STD of amplitude The upper peak amplitude, mean amplitudeand STD of amplitude of the upper peaks in predefined time window. Theupper peak represent the depth of respiration. 3 35 Respi- ratory LowerPeak amplitude, mean amplitude and STD of amplitude The lower peakamplitude, mean amplitude and STD of amplitude of the lower peaks inpredefined time window. The lower peaks represent the depth of breathrelease. 2 36 Respi- ratory Respiratory rate, mean rate and std rate Therate, mean rate and STD of rate in predefined time window. The rate is1/Peak-to-peak interval. 3 37 BP BP Peak (P) and Trough (T) amplitude,mean amplitude and std of amplitude The amplitude of the Peak (P) andthe Trough (T) of the NIBP signal, mean amplitudes and STD of amplitudesin predefined time window. Peak denotes a systolic blood pressure;Trough denotes a diastolic blood pressure. 6 38 EEG/ EMG Power of A, β,γ, δ, θ frequency bands Power of the frequency bands of power spectrumof EEG/fEMG signal in predefined time segment. Frequency range ofdifferent bands [Hz]: delta, δ 0.5-4 theta, θ 4-8 alpha, α 8-14 beta, β14-30 gamma γ, 30-70 4 39 EEG/ EMG Mean frequency The sum of the productof the power spectrum values in predefined time segment and thefrequencies, divided by the total power. 1 40 EEG/ EMG Peak frequencyThe frequency of the highest peak of the power spectrum in predefinedtime segment. 1 41 EEG/ EMG Spectral Edge Frequency The frequency belowwhich x percent of the power are located. Typically x is in the range 75to 95. 1 42 EEG/ EMG Approximate Entropy- For details see (Bruhn, Ropckeand Hoeft 2000) 1 43 EEG/ EMG BSR-Burst Suppression ratio$\begin{matrix}{{The}\mspace{14mu}{burst}\mspace{14mu}{suppression}\mspace{14mu}{ratio}} \\{{is}\mspace{14mu}{the}\mspace{14mu}{proportion}\mspace{14mu}{of}\mspace{14mu}{the}} \\{{suppression}\mspace{14mu}{EEG}\mspace{14mu}{in}\mspace{14mu}{the}} \\{{analyzed}\mspace{14mu}{epoch}} \\{\left( {{usually}\mspace{14mu}{one}\mspace{14mu}{minute}} \right)\text{:}} \\{{BSR} = {\frac{\begin{matrix}{{total}\mspace{14mu}{time}\mspace{14mu}{of}} \\{suppression}\end{matrix}}{{epoch}\mspace{14mu}{length}}100}}\end{matrix}\quad$ 1 44 EEG/ EMG BcSEF $\begin{matrix}{{Burst}\mspace{14mu}{compensated}\mspace{14mu}{spectral}} \\{{edge}\mspace{14mu}{frequency}} \\{{BcSEF} = {{SEF}\mspace{11mu}\left( {1 - \frac{{BSR}\mspace{11mu}\%}{100\%}} \right)}}\end{matrix}\quad$ 1 45 EEG/ EMG WSMF A generalized form of spectraledge frequency, referred to as weighted spectral median frequency(WSMF), edge frequency is calculated not necessarily from PSD but fromamplitude spectrum, which is raised to the power p = [0.1 . . . 2.4];second, the cutoff frequencies of the original spectrum arewell-defined; and, third, factor r = [0:05 : : : 0:95] is used, thepercentile of the spectrum (e.g., r = 0:5 for MF and r = 0:95 for SEF).1 46 EEG/ EMG CUP $\begin{matrix}{{Canonical}\mspace{14mu}{univariate}} \\{{parameter}\text{:}\mspace{14mu}{frequency}\mspace{14mu}{bins}} \\{{with}\mspace{14mu} a\mspace{14mu}{width}\mspace{14mu}{of}\mspace{14mu} 3\mspace{14mu}{Hz}\mspace{14mu}{or}} \\{{classical}\mspace{14mu}{frequency}\mspace{14mu}{bands}} \\{{are}\mspace{14mu}{optimally}\mspace{14mu}{weighted}\mspace{14mu}{to}} \\{{obtain}\mspace{14mu}{the}\mspace{14mu}{best}\mspace{14mu}{possible}} \\{{correlation}\mspace{14mu}{with}\mspace{14mu}{the}\mspace{14mu}{{drugs}'}} \\{{effect}\text{-}{site}\mspace{14mu}{concentration}\mspace{14mu}{as}} \\{{{obtained}\mspace{14mu}{from}}\mspace{14mu}} \\{{pharmacokinetic}\text{-}} \\{{pharmacodynamic}\mspace{14mu}\left( {{PK}\text{-}{PD}} \right)} \\{modeling} \\{{CUP} = {\sum\limits_{k = 1}^{10}{\gamma_{k}\log\mspace{11mu} p_{k}}}}\end{matrix}\quad$ 1 47 EEG/ EMG SpEn- $\begin{matrix}{{Spectral}\mspace{14mu}{Entropy}} \\{{SpEn} = {- {\sum\limits_{k}^{N}{p_{k}\log\mspace{11mu}{p_{k}.}}}}}\end{matrix}\quad$ 1 48 EEG/ EMG BcSpEn- $\quad\begin{matrix}{{Burst}\mspace{14mu}{compensated}\mspace{14mu}{Spectral}} \\{Entropy} \\{{BcSpEn} = {{{SpEn}\begin{pmatrix}{1 -} \\\frac{{BSR}\mspace{11mu}(\%)}{100\%}\end{pmatrix}}.}}\end{matrix}$ 1 49 EEG/ EMG Beta Ratio${BetaRatio} = {\log\;{\frac{{\hat{P}}_{30\text{-}47{Hz}}}{{\hat{P}}_{11\text{-}20{Hz}}}.}}$1 50 EEG/ EMG Histogram parameters Mean, Standard deviation, Kurtosis,Skewness of signal histogram in predefined time segment. 4 51 EEG/ EMGAR parameters Parameters of AR representation (Schlogl 2006) N 52 EEG/EMG Normalized slope descriptors (Hjorth parameters) NSD parameters canbe defined by means of first and second derivatives. “Activity” is ameasure of the mean power, “Mobility” is an estimate of the meanfrequency and “Complexity” is an estimate of the bandwidth of the signal(frequency spread) (Hjorth 1973). 3 53 EEG/ EMG Barlow parametersParameters based on Barlow EEG model which is an alternative timefrequency decomposition. Parameters such as Running Mean Frequency andSpectral Purity Index (Goncharova and Barlow 1990) 3 54 EEG/ EMGWackermann parameters Three multi-channel linear descriptors of EEGsignal. spatial complexity (Ω), field power (Σ) and frequency of fieldchanges (Φ) (Wackermann 1999) 3 55 EEG/ EMG Brain rate Weighted MeanFrequency (Pop-Jordanova and Pop- Jordanov 2005) 1 56 EEG/ EMGSynchFastSlow $\begin{matrix}{{SynchFastSlow} = {\log\;{\frac{{\hat{B}}_{40\text{-}47{Hz}}}{{\hat{B}}_{0.5\text{-}47{Hz}}}.}}} \\{{{The}\mspace{14mu}{spectrum}\mspace{14mu}{and}\mspace{14mu}{bispectrum}},} \\{{derived}\mspace{14mu}{from}\mspace{14mu}{two}\text{-}{second}} \\{{epochs},{{are}\mspace{14mu}{smoothed}\mspace{14mu}{using}\mspace{14mu} a}} \\{{running}\mspace{14mu}{average}\mspace{14mu}{against}\mspace{14mu}{those}} \\{{calculated}\mspace{14mu}{in}\mspace{14mu}{the}\mspace{14mu}{previous}} \\{{{minute}.\mspace{14mu} 3}\mspace{14mu}{minutes}\mspace{14mu}{window}} \\{{is}\mspace{14mu}{required}\mspace{14mu}{to}\mspace{14mu}{obtain}\mspace{14mu} a} \\{{consistent}\mspace{14mu}{estimate}\mspace{14mu}{of}\mspace{14mu}{the}} \\{{bicoherence}.}\end{matrix}\quad$ 1 57 EEG/ EOG 80 Hz frequency in EEG near the eyesOcular microtremor (OMT) is a constant, physiological, high frequency(peak 80 Hz), low amplitude (estimated circa 150- 2500 nm) eye tremor. 158 EMG Spectrum analysis- Power of the frequency bands of power spectrumof EMG signal in predefined time segment. Frequency range of differentbands [Hz]: 1 59 EMG mean frequency The sum of the product of the powerspectrum values in predefined time segment and the frequencies, dividedby the total power. 1 60 EMG Peak frequency The frequency of the highestpeak of the power spectrum in predefined time segment. 1 61 EMG Totalpower The sum of the power spectrum within the epoch 1 62 EMGSpontaneous lower oesophageal contractions (SLOC) Lower oesophagealcontractility (LOC). Spontaneous lower oesophageal contractions (SLOC)are non-propulsive spontaneous contractions mediated via vagal motornuclei and reticular activating system in the brain stem. The frequencyof these movements is increased as the dose of the anaesthetic isreduced. (Thomas and Evans 1989) 1 63 Airway CO2 Average/ VariabilityEnd tidal Carbon Dioxid (anesthesia) 2 64 Airway Gases Average End tidalsevofluane (anesthesia) 1 65 acceler- ometer X, Y, Z θ Average value,Variability accelerometer X, Y, Z theta, movement analysis 12

Within the context of the present invention the term normalizationrefers to a signal processing or preprocessing technique at is known andaccepted in the art for manipulating a physiological signal and/orfeature. Normalization process for example comprises but is not limitedto feature normalization to removes a subject's baseline; featurenormalization to normalize a subject's baseline feature variability,feature normalization to remove a subject's baseline feature mean,feature normalization which removes the feature mean and normalize thefeature variability; feature normalization which normalize the featurevalue into a value between [0,1]. Most preferably normalization isperformed following feature extraction and is proceeded by at least oneor more signal preprocessing stages for example including but notlimited to raw data preprocessing, optionally by noise filtering and/orartifact reduction, then features or the like are extracted and finallynormalized. Optionally and preferably, normalization may be utilized foreither a historical data set or the current data set. Optionallynormalization may be performed on any data, parameters, variables and orfeatures associated with the method and system of the present inventionfor example including but not limited to kurtosis, skewness, higherorder moments and cummulants, probability distribution functionsassociated with a feature or the like.

Within the context of the present invention the term feature selectionand dimensionality reduction refers to the process and/or technique ofidentifying, abstracting, representing or manipulating a plurality ofphysiological features and/or physiological signals in a more compactand/or reduced form. Optionally, in the process of feature selection anddimensionality reduction the number of physiological features and/orsignals used with the system and method of the present invention isreduced.

Optionally and preferably feature selection and dimensional reductiontechniques optionally employed within the system and method of thepresent invention may be linear or nonlinear, for example including butare not limited to Multi Dimensional Scaling (Borg, et. al., 2005),Principal Component Analysis (PCA), Sparse PCA (SPCA), Fisher LinearDiscriminant Analysis (FLDA), Sparse FLDA (SFLDA), Kernel PCA (KPCA)(Scholkopf, et al., 1998), ISOMAP (Tenenbaum, et al., 2000), LocallyLinear Embedding (LLE) (Roweis, et al., 2000), Laplacian Eigenmaps(Belkin, et al., 2003), Diffusion Maps (Coifman, et al., 2005), HessianEigenmaps (Donoho, et al., 2003), Independent Component Analysis (ICA),Factor analysis (FA), Hierarchical Dimensionality Reduction (HDR), SureIndependence Screening (SIS), Fisher score ranks, t-test rank,Mann-Whitney U-test taken alone or in any combination thereof or thelike feature selection and dimensionality reduction techniques as isknown and accepted in the art.

Within the context of the present invention the term classificationrefers to the analysis performed according to the system and method ofthe present invention to identify and or determine the physiologicalstate of a subject with respect to pain. Optionally, the classificationmay be rendered in at least two and optionally a plurality of classesassociated with pain. Optionally and preferably a plurality ofclassification techniques also optionally referred to herewith asclassifiers may be utilized according to the system and method of thepresent invention as is known and accepted in the art or in anycombination thereof for example including but is not limited to NearestShrunken Centroids (NSC) (Tibshirani, et al., 2002), Classification andRegression Trees (CART) (Hastie, et al., 2009), ID3, C4.5, MultivariateAdditive regression splines (MARS), Multiple additive regression trees(MART), Nearest Centroid (NC) classifier (Hastie, et al., 2009),Shrunken Centroid Regularized Linear Discriminate and Analysis (SCRLDA)(Guo, et al., 2007), Random Forest (Breiman, 2001), Boosting (Hastie, etal., 2009), Bagging Classifier (Breiman, 1996), AdaBoost, RealAdaBoost,LPBoost, TotalBoost, BrownBoost, MadaBoost, LogitBoost, GentleBoost,RobustBoost, Support Vector Machine (SVM) (Vapnik, 1998), kernelizedSVM, Linear classifier, Quadratic Discriminant Analysis (QDA)classifier, Naïve Bayes Classifier and Generalized Likelihood Ratio Test(GLRT) classifier with plug-in parametric or non-parametric classconditional density estimation, k-nearest neighbor, Radial Base Function(RBF) classifier, Multilayer Perceptron classifier, Bayesian Network(BN) classifier (Hastie, et al., 2009) (Bishop, 2006) or the like as areknown and accepted in the art.

Optionally and preferably classification is performed using a classifieradept at multi-class classification. Optionally, multi-classclassification may be adapted from binary classifier as is known andaccepted in the art. Optionally binary classifiers may be adapted toperform a multi-class classification by reducing the multi-class problemto a plurality of multiple binary problems using methods as is known andaccepted in the art for example including but not limited to one-vs-onewith voting schemes by majority vote or pairwise coupling, one-vs-rest,Error Correcting Output Codes, or the like as is described and known inthe art (Bishop, 2006) (Hastie, et al., 2009).

Optionally hierarchical multi-class classification may be performed as atree structure having a single parent class or directed acyclic graphstructure having at least one or more parental class.

Most preferably, the classifier is trained on a “training set”.Optionally the training set may comprise data for classification that ismade available from a plurality of sources, for example including butnot limited to publicly available databases, proprietary clinical trialsdata, on site recorded data from at least one or more subject. Mostpreferably, the training set comprises input and output signals thatmimic the input and output signals of the pain classifier according tothe present invention. Most preferably the input signals comprising thetraining set are similar in nature to the expected input, of a painclassifier according to the present invention. Optionally and preferablythe training set input signals comprise data similar to the featuresassociated with a extracted feature vector according to the presentinvention, and/or a priori data as described above. Most preferably theoutput signals comprising the training set are similar in nature to theexpected output from a pain classifier, according to the presentinvention. Most preferably, the training set is compiled by a painexpert for example a physician or other skilled to person in the art ofpain detection. Optionally, the training set is complied during aclinical trial comprising controlled pain stimuli.

Within the context of the present invention the term learning and ortraining refers to the process of training a pain classifier based on agiven training, as is known and accepted in the art. Most preferably,the training process will provide for classifying previously unseeninput data, not from the training set, with sensitivity and specificityoptionally similar and more preferably better than the performances of ahuman operator.

A preferred embodiment of the present invention provides for a systemand method for detecting, classifying and monitoring pain. Mostpreferably, the system and method of the present invention utilizes acollection of features referred to as a great plurality of features,hereinafter referred to as GPF, optionally comprising at least 9physiological features, for detecting, classifying and identifying painlevel. Optionally and preferably, a priori data, as described above, maybe coupled with the GPF to facilitate pain monitoring and detection.

Optionally a plurality of features (GPF) providing for painclassification according to embodiments of the present invention maycomprise a plurality of features corresponding to and extracted from atleast three or more physiological signals, as depicted in Table 1.

Preferably a first vector comprising a plurality of features (GPF)includes features, and/or feature groups originating from at least threephysiological signals including PPG, GSR and at least one or morephysiological signals chosen from the group consisting: of skintemperature, ECG, Respiration, EMG, and EEG/FEMG. Optionally and mostpreferably a first vector of features may be built from the featuresassociated with at least 1 (one) of 32 (thirty two) optionalphysiological signal groups, for example including but not limited to:

PPG-GSR; PPG-GSR-TEMP; PPG-GSR-ECG; PPG-GSR-RESP; PPG-GSR-EMG;PPG-GSR-EEG/FEMG; PPG-GSR-TEMP-ECG; PPG-GSR-TEMP-RESP; PPG-GSR-TEMP-EMG;PPG-GSR-TEMP-EEG/FEMG; PPG-GSR-ECG-RESP; PPG-GSR-ECG-EMG;PPG-GSR-ECG-EEG/FEMG; PPG-GSR-RESP-EMG; PPG-GSR-RESP-EEG/FEMG;PPG-GSR-EMG-EEG/FEMG; PPG-GSR-TEMP-ECG-RESP; PPG-GSR-TEMP-ECG-EMG;PPG-GSR-TEMP-ECG-EEG/FEMG; PPG-GSR-TEMP-RESP-EMG;PPG-GSR-TEMP-RESP-EEG/FEMG; PPG-GSR-ECG-RESP-EMG;PPG-GSR-ECG-RESP-EEG/FEMG; PPG-GSR-ECG-EMG-EEG/FEMG;PPG-GSR-RESP-EMG-EEG/FEMG; PPG-GSR-TEMP-EMG-EEG/FEMG;PPG-GSR-TEMP-ECG-RESP-EMG; PPG-GSR-TEMP-ECG-RESP-EEG/FEMG;PPG-GSR-TEMP-ECG-EMG-EEG/FEMG; PPG-GSR-TEMP-RESP-EMG-EEG/FEMG;PPG-GSR-ECG-RESP-EMG-EEG/FEMG; PPG-GSR-TEMP-ECG-RESP-EMG-EEG/FEMG; orthe like group or combination of physiological signals.

Optionally and preferably each of the physiological signals namely, PPG,GSR, skin temperature, ECG, Respiration, EMG, and EEG/FEMG, or the likemay individually provide a feature and/or a group of features based onthe physiological signal itself, each feature group defined as followsaccording to its physiological signal, as follows:

PPG features, for example include but not limited to:

PPG Peak (P) amplitude, Trough (T) amplitude, mean PPG Peak (P)amplitude, and std of PPG Peak (P) amplitude, mean Trough (T) amplitude,and std of Trough (T) amplitude, PPG peak to peak time intervals, PPGpeak to peak interval mean and PPG peak to peak interval std; powerspectrum of the PPG peak to peak intervals: VLF Power, LF Power and HFPower;

GSR features, table 1, for example including but not limited to:

GSR amplitude, GSR mean amplitude and GSR amplitude std, GSR Peak (P)amplitude, mean Peak (P) amplitude and Peak (P) amplitude std; GSR peakto peak time intervals, mean GSR peak to peak time interval; and GSRpeak to peak time intervals std; Phasic EDA: amplitude, mean amplitudeand std of amplitude;

Skin temperature features, Table 1, for example including but notlimited to:

Temperature amplitude, mean amplitude and std of amplitude; Temp Peak(P) amplitude, mean amplitude and std of amplitude; Temperature peak topeak time intervals, mean and std (variability) of interval;

ECG features, Table 1, for example including but not limited to:

ECG-PPG PFT Pulse Transition time; ECG R to R time intervals, mean andstd (variability) of intervals; Power of VLF, LF and HF frequency bandsof power spectrum of the ECG R to R intervals (heart rate variability);

Respiration features, Table 1, for example including but not limited to:Upper peak amplitude, mean amplitude and STD of amplitude, Respiratoryrate, mean rate and std rate

EMG features for example including but not limited to: Power of thefrequency bands of power spectrum of EMG signal; EMG Power Spectrum Meanfrequency;

EMG Power Spectrum Highest Peak Frequency;

EEG/FEMG features, table 1, for example including but not limited to:Power of the alpha (α), beta (β), gamma (γ), delta (δ), theta (θ)frequency bands of power spectrum of EEG/FEMG signal; EMG Power SpectrumMean frequency; EMG Power Spectral edge frequency ‘Coherence between 2or more EEG/FEMG channels;

Preferably any combination of the feature groups, defined above, andmost preferably a combination of features according to at least one ormore of the above referenced 32 (thirty two) physiological signal groupsmay be utilized to form a first vector of feature upon whichclassification will be based and performed according to embodiments ofthe present invention.

Optionally the system and method of the present invention provide forclassification of pain into a plurality of classes, most preferablycomprising at least 2 classes. Optionally, pain may be categorized into3 or more classification groups. Optionally, each pain classificationmay be further provided with scalable scoring method optionally using anumerical scoring for example including a scale of 1 to 10, a scale of1-100 or the like. Optionally, pain classification scoring may befurther correlated to a subjective pain scaling and/or scoring schemesfor example including but not limited to Visual Analog Scale (VAS),Numeric Pain Scale (NPS), verbal scale or the like as is known in theart. Optionally, pain may be categorized according to a disease,stimulus and or medicament for example pain associated with cancer wouldbe classified in one class while pain associated with diabetes would beclassified in another class.

Optionally, the system and method according to the present invention maybe adapted to and/or provide for pain monitoring, classification and/ordetection according to at least one or more optional sources of pain forexample including but not limited to stimulus, medicament or a disease.

An optional embodiment of the system and method according to the presentinvention is provided with a GPF that are most preferably abstractedfrom a plurality of physiological signals for example including but notlimited to blood to pressure, respiration, internal or skin temperature,pupil diameter, GSR, and signals received from ECG, PPG, EOG, EGG, EEG,EMG, EGG, LDV, capnograph and accelerometer or the like physiologicalsignal as is known and accepted in the art.

Optionally and preferably, GPF is abstracted from a plurality ofphysiological signals using at least one or more feature extractiontechniques as described above and as is known and accepted in the art.

Optionally, the GPF according to the present invention may be furtherprocessed using at least one or more feature selection anddimensionality reducing techniques as described above and as is knownand accepted in the art, most preferably providing a set of parametersbased on the GPF that are preferably then utilized for classification.

Most preferably, classification is provided by classifiers as describedabove and as is known and accepted in the art optionally, the classifiermay provide for linear and/or non-linear classification. Optionally theclassifier provides for classification of at least 2 classes. Optionallyand preferably the classifier provides for multi-class classification ofa plurality of classes

An optional embodiment for a system for pain monitoring and detectionaccording to the present invention comprises a signal acquisitionmodule, a processing module and a communication module. Optionally thesystem may further comprise a display module. Most preferably, signalacquisition comprises a plurality of optional sensors and/or transducersfor measuring and or obtaining physiological signals and data. Mostpreferably, processing module provides for processing the physiologicalsignals provided through the acquisition module, for example to abstractfrom the signals the GPF. Optionally, the processing module provide forall processing of signals for example including but not limited tofeature abstraction, preprocessing, feature selection and dimensionalityreduction and classification. Optionally and preferably, thecommunication module is provided to communicate the classificationresults from the processing module. Optionally the communication modulemay provide for communicating results to an auxiliary person, system,device, machine, processor for example including but not limited to ahigher processing center, person, caregiver, call center, or the like inany combination thereof.

Optionally, the display module is provided to display and/or communicatethe classification results from the processing module. Optionally,display module comprise a display for example including but not limitedto a visual display, printed display or audible display, or the like fordisplaying and conveying pain monitoring. Optionally, processing modulemay be realized as a wired or wireless device for example including butnot limited to a computer, a server, PDA, mobile telephone, display,printout or the like as is known and accepted in the art.

An optional embodiment according to the present invention provides for amethod for detecting and classifying the pain status of a patient byanalyzing a plurality of physiological signals, the method comprises:

a. Signal acquisition for acquiring the plurality of physiologicalsignals; and

b. pre-processing the acquired plurality of physiological signals toimprove signal quality comprising at least one or more chosen from thegroup consisting of synchronization, noise filtering, artifactreduction, therein forming a plurality of pre-processed physiologicalsignals; and

c. processing the pre-processed plurality of physiological signals, theprocessing comprising:

-   -   i. feature extraction from at least three or more physiological        signals including PPG, GSR and at least one or more        physiological signals chosen from the group consisting: of skin        temperature, ECG, Respiration, EMG, and EEG/FEMG to facilitate        detection of pain in a patient; and forming a first vector        comprising a set of extracted features; and    -   ii. transforming the first vector to a second vector wherein        pain detection is performed based on the second vector and        wherein the transformation comprises normalization, feature        selection and dimensionality reduction techniques; and    -   iii. detecting the pain status of a patient by applying a        classification function to classify the second vector into at        least two classes of pain.

Optionally the method further comprises communicating the detected painstatus of the patient to at least one or more for example including butnot limited to a higher processing center, person, caregiver, callcenter and any combination thereof.

Optionally and preferably the first vector comprises a plurality offeatures including parameters extracted from the PPG and GSRphysiological signals including:

a. PPG features chosen from the group consisting of:

-   -   i. PPG Peak (P) amplitude, Trough (T) amplitude, mean PPG        Peak (P) amplitude, and std of PPG Peak (P) amplitude, mean        Trough (T) amplitude, and std of Trough (T) amplitude; and    -   ii. PPG peak to peak time intervals, PPG peak to peak interval        mean and PPG peak to peak interval std; and    -   iii. power spectrum of the PPG peak to peak intervals: VLF        Power, LF Power and HF Power; and

b. GSR features chosen from the group consisting of:

-   -   i. GSR amplitude, GSR mean amplitude and GSR amplitude std; and    -   ii. GSR Peak (P) amplitude, mean Peak (P) amplitude and Peak (P)        amplitude std; and    -   iii. GSR peak to peak time intervals, mean GSR peak to peak time        interval; and GSR peak to peak time intervals std; and    -   iv. Phasic EDA: amplitude, mean amplitude and std of amplitude.

Optionally, the first vector further comprises a group of featuresextracted from one of the physiological signals selected from the groupconsisting of skin temperature, ECG, Respiration, EMG, and EEG/FEMG.

Optionally, the first vector further comprises a group of featuresextracted from two of the physiological signals selected from the groupconsisting of skin temperature, ECG, Respiration, EMG, and EEG/FEMG.

Optionally, the first vector further comprises a group of featuresextracted from three of the physiological signals selected from thegroup consisting of skin temperature, ECG, Respiration, EMG, andEEG/FEMG.

Optionally, the first vector further comprises a group of featuresextracted from four of the physiological signals selected from the groupconsisting of skin temperature, ECG, Respiration, EMG, and EEG/FEMG.

Optionally, the first vector further comprises a group of featuresextracted from the physiological signals essentially consisting of PPG,GSR, skin temperature, ECG, Respiration, EMG, and EEG/FEMG.

Optionally and preferably, extracted features are chosen form at leastone or more of the feature groups corresponding to a physiologicalsignal chosen from the group consisting of skin temperature, ECG,Respiration, EMG, and EEG/FEMG, the feature groups comprising:

a. Skin temperature features chosen from the group of consisting:

-   -   i. Temperature amplitude, mean amplitude and std of amplitude;        and    -   ii. Temp Peak (P) amplitude, mean amplitude and std of        amplitude; and    -   iii. Temperature peak to peak time intervals, mean and std        (variability) of interval; and

b. ECG features chosen from the group of consisting:

-   -   i. ECG-PPG PTT Pulse Transition time; and    -   ii. ECG R to R time intervals, mean and std (variability) of        intervals; and    -   iii. Power of VLF, LF and HF frequency bands of power spectrum        of the ECG R to R intervals (heart rate variability); and

c. Respiration features chosen from the group of consisting of:

-   -   i. Upper peak amplitude, mean amplitude and STD of amplitude;        and    -   ii. Respiratory rate, mean rate and std rate; and

d. EMG features chosen from the group of consisting of:

-   -   i. Power of the frequency bands of power spectrum of EMG signal;        and    -   ii. EMG Power Spectrum Mean frequency; and    -   iii. EMG Power Spectrum Highest Peak Frequency; and

e. EEG/FEMG features chosen from the group of consisting of:

-   -   i. Power of the alpha (α), beta (β), gamma (γ), delta (δ), theta        (θ) frequency bands of power spectrum of EEG/FEMG signal; and    -   ii. EMG Power Spectrum Mean frequency; and    -   iii. EMG Power Spectral edge frequency; and    -   iv. Coherence between 2 or more EEG/FEMG channels.

Optionally, signal acquisition facilitates acquiring physiologicalsignals with a plurality of sensors for example including but notlimited to: ECG, PPG, blood pressure, respiration, internal bodytemperature, skin temperature, EOG, pupil diameter monitoring, GSR,EEG/FEMG, EMG, EGG, LDV, capnograph and accelerometer, or the like.Optionally and preferably, the physiological signals are pre-processedand processed, facilitating said detection and classification of painstatus.

Optionally the method further comprises, obtaining and processing apriori data for facilitating the detection and classification of painstatus. Optionally a priori data may for example include but is notlimited to patient associated data, historical data, data from the groupconsisting of data supplied by the physician's, environmentalparameters, patient parameters, disease, stimulus and medicament, anycombination thereof, or the like.

Optionally, the second vector further comprises a priori data.

Optionally, pain classification may be classified into groups forexample including but not limited to: at least two classes, at leastthree classes, pain/no-pain states, graduated scale, scale of 1 to 10, ascale of 1-100, subjective pain scaling schemes, subjective scoringschemes, Visual Analog Scale (VAS), Numeric Pain Scale (NPS) and verbalscale, or the like.

Optionally, dimensionality reducing transformation comprises a linear ornon-linear transformation, for example including but is not limited to:Multi Dimensional Scaling (MDS), Principal Component Analysis (PCA),Sparse PCA (SPCA), Fisher Linear Discriminant Analysis (FLDA), SparseFLDA (SFLDA), Kernel PCA (KPCA), ISOMAP, Locally Linear Embedding (LLE),Laplacian Eigenmaps, Diffusion Maps, Hessian Eigenmaps, IndependentComponent Analysis (ICA), Factor analysis (FA), Hierarchical Featureselection and dimensionality reduction (HDR), Sure IndependenceScreening (SIS), Fisher score ranks, t-test rank, Mann-Whitney U-test orany combination thereof.

Optionally classification function may be a linear or non-lineartransformation for example including but not limited to: NearestShrunken Centroids (NSC), Classification and Regression Trees (CART),ID3, C4.5, Multivariate Additive regression splines (MARS), Multipleadditive regression trees (MART), Nearest Centroid (NC), ShrunkenCentroid Regularized Linear Discriminate and Analysis (SCRLDA), RandomForest, Boosting, Bagging Classifier, AdaBoost, RealAdaBoost, LPBoost,TotalBoost, BrownBoost, MadaBoost, LogitBoost, GentleBoost, RobustBoost,Support Vector Machine (SVM), kernelized SVM, Linear classifier,Quadratic Discriminant Analysis (QDA) classifier, Naïve Bayes Classifierand Generalized Likelihood Ratio Test (GLRT) classifier with plug-inparametric or non-parametric class conditional density estimation,k-nearest neighbor, Radial Base Function (RBF) classifier, MultilayerPerceptron classifier, Bayesian Network (BN) classifier, multi-classclassifier adapted from binary classifier with one-vs-one majorityvoting, one-vs-rest, Error Correcting Output Codes, hierarchicalmulti-class classification, Committee of classifiers, or the like, andany combination thereof.

Optionally, the method according to an optional embodiment of thepresent invention comprises: feature selection and dimensionalityreduction provided by Hierarchical Dimensionality Reduction (HDR); andclassification provided by a Random Forest classifier.

Optionally, the method according to an optional embodiment of thepresent invention comprises, feature selection and dimensionalityreduction provided by Fisher Score rank, and SFLDA and classificationprovided by a RealAdaboost classifier within a Boosting framework.

Optionally classification may be adapted for pain experienced with aparticular disease, stimulus or medicament.

Optionally, the status of pain may be adapted for a patient in variousstates of consciences chosen from the group consisting of: sedated,partially sedate and awake, semi-awake.

An optional embodiment according to the present invention provides for asystem for detecting and classifying the pain status of a patient byanalyzing a plurality of physiological signals comprising:

a. a signal acquisition module comprising a plurality of sensors and/ortransducers for measuring and/or obtaining the at least three or morephysiological signals and a priori data from a subject; and

b. a processing module for processing the physiological signals,comprising:

-   -   i. pre-processing the acquired plurality of physiological        signals to improve signal quality comprising at least one or        more tools for example including but not limited to        synchronization, noise filtering, artifact reduction therein        forming a pre-processed plurality of physiological signals; and    -   ii. processing the pre-processed plurality of physiological        signals comprising:    -   iii. feature extraction from at least three or more        physiological signals including PPG, GSR and at least one or        more physiological signals chosen from the group consisting: of        skin temperature, ECG, Respiration, EMG, and EEG/FEMG to        facilitate detection of pain in a patient; and forming a first        vector comprising a set of extracted features; and    -   iv. transforming the first vector to a second vector wherein        pain detection is performed based on the second vector and        wherein the transformation comprises normalization and feature        selection and dimensionality reduction techniques; and    -   v. detecting the pain status of a patient by applying a        classification function to classifying the second vector into at        least two classes of pain; and

c. a communicating module for communicating the detected pain status ofthe patient to at least one or more for example including but notlimited to a higher processing center, person, caregiver, call centerand any combination thereof.

Optionally, the system may further comprise a display module fordisplaying the detected and classified pain stimulus.

Optionally and preferably, the processing module provides for extractingfeatures form at least one or more of the feature groups correspondingto a physiological signal chosen from the group consisting of skintemperature, ECG, Respiration, EMG, and EEG/FEMG, the feature groupscomprising:

a. Skin temperature features chosen from the group of consisting:

-   -   i. Temperature amplitude, mean amplitude and std of amplitude;        and    -   ii. Temp Peak (P) amplitude, mean amplitude and std of        amplitude; and    -   iii. Temperature peak to peak time intervals, mean and std        (variability) of interval; and

b. ECG features chosen from the group of consisting:

-   -   i. ECG-PPG PTT Pulse Transition time; and    -   ii. ECG R to R time intervals, mean and std (variability) of        intervals; and    -   iii. Power of VLF, LF and HF frequency bands of power spectrum        of the ECG R to R intervals (heart rate variability); and

c. Respiration features chosen from the group of consisting of:

-   -   i. Upper peak amplitude, mean amplitude and STD of amplitude;        and    -   ii. Respiratory rate, mean rate and std rate; and

d. EMG features chosen from the group of consisting of:

-   -   i. Power of the frequency bands of power spectrum of EMG signal;        and    -   ii. EMG Power Spectrum Mean frequency; and    -   iii. EMG Power Spectrum Highest Peak Frequency; and

e. EEG/FEMG features chosen from the group of consisting of:

-   -   i. Power of the alpha (α), beta (β), gamma (γ), delta (δ), theta        (θ) frequency bands of power spectrum of EEG/FEMG signal; and    -   ii. EMG Power Spectrum Mean frequency; and    -   iii. EMG Power Spectral edge frequency; and    -   iv. Coherence between 2 or more EEG/FEMG channels.

Optionally the system according to the present invention furthercomprises a signal acquisition module that facilitates acquiringphysiological signals with a plurality of sensors for example includingbut not limited to: ECG, PPG, blood pressure, respiration, internal bodytemperature, skin temperature, EOG, pupil diameter monitoring, GSR,EEG/FEMG, EMG, EGG, LDV, capnograph and accelerometer, or the like.Optionally, the physiological signals may be pre-processed and processedfacilitating the detection and classification of pain status.

Unless otherwise defined the various embodiment of the present inventionmay be provided to an end user in a plurality of formats, platforms, andmay be outputted to at least one of a computer readable memory, acomputer display device, a printout, a computer on a network or a user.

Unless otherwise defined, all technical and scientific terms used hereinhave the same meaning as commonly understood by one of ordinary skill inthe art to which this invention belongs. The materials, methods, andexamples provided herein are illustrative only and not intended to belimiting.

Implementation of the method and system of the present inventioninvolves performing or completing certain selected tasks or stepsmanually, automatically, or a combination thereof. Moreover, accordingto actual instrumentation and equipment of preferred embodiments of themethod and system of the present invention, several selected steps couldbe implemented by hardware or by software on any operating system of anyfirmware or a combination thereof. For example, as hardware, selectedsteps of the invention could be implemented as a chip or a circuit. Assoftware, selected steps of the invention could be implemented as aplurality of software instructions being executed by a computer usingany suitable operating system. In any case, selected steps of the methodand system of the invention could be described as being performed by adata processor, such as a computing platform for executing a pluralityof instructions.

Although the present invention is described with regard to a “computer”on a “computer network”, it should be noted that optionally any devicefeaturing a data processor and/or the ability to execute one or moreinstructions may be described as a computer, including but not limitedto a PC (personal computer), a server, a minicomputer, a cellulartelephone, a smart phone, a PDA (personal data assistant), a pager. Anytwo or more of such devices in communication with each other, and/or anycomputer in communication with any other computer, may optionallycomprise a “computer network”.

BRIEF DESCRIPTION OF THE DRAWINGS

The invention is herein described, by way of example only, withreference to the accompanying drawings. With specific reference now tothe drawings in detail, it is stressed that the particulars shown are byway of example and for purposes of illustrative discussion of thepreferred embodiments of the present invention only, and are presentedin order to provide what is believed to be the most useful and readilyunderstood description of the principles and conceptual aspects of theinvention. In this regard, no attempt is made to show structural detailsof the invention in more detail than is necessary for a fundamentalunderstanding of the invention, the description taken with the drawingsmaking apparent to those skilled in the art how the several forms of theinvention may be embodied in practice.

In the drawings:

FIG. 1 is schematic block diagram of an exemplary system according tothe present invention.

FIGS. 2A-C are schematic block diagram of the machine learning module ofFIG. 1 in greater detail.

FIGS. 3A-E are schematic block diagrams of optional embodiments of asystem according to the present invention.

FIG. 4 is an exemplary method according to the present invention.

FIGS. 5A-C are scatter plots depicting pain classification according toan optional system and method of the present invention.

DESCRIPTION OF THE PREFERRED EMBODIMENTS

The present invention is of a system and a method for pain detection andmonitoring most preferably detection and monitoring of pain isfacilitated by processing physiological signals and features mostpreferably represented by a data vector comprising a great plurality offeatures. The principles and operation of the present invention may bebetter understood with reference to the drawings and the accompanyingdescription.

For the sake of clearly throughout the figures similar labels andnumbering scheme is used throughout for equivalent or similarlyfunctioning elements.

FIG. 1 is a schematic block diagram of an exemplary system 100 accordingto the present invention for pain monitoring comprising physiologicalsignal acquisition module 102 and processing module 104.

Most preferably signal acquisition module 102 is provided for acquiringand measuring physiological measurements from a subject able toexperience pain optionally including but not limited to a person isoptionally provided with at least one or more preferably a plurality ofsensor or transducers as are known and accepted in the art. Processingmodule 104 is most preferably provided with a plurality of sub-modulesfor example including but not limited to signal processing module 2,feature extraction module 3, machine learning module 4. Optionally andpreferably signal processing module 2 provides for preprocessing of theacquired signal for example including but not limited to normalization,filtering, noise reduction, SNR optimization, domain transformations,statistical analysis, spectral analysis, wavelet analysis, or the like.

Optionally and preferably feature extraction module 3 is provided forextracting features associated with the signals acquired with module102. Optionally and preferably, feature extraction module may be furtherprovided with an a priori data sub-module 6 most preferably forproviding data beyond the acquired signals for example including but notlimited to physician and/or caregiver data, medical history, geneticpredispositions data, medical records or the like a priori data. Mostpreferably feature extraction is provided by performing signalprocessing techniques as is known and accepted in the art to extractfrom a physiological signal relevant and pertinent data. For example, aECG measured acquired in module 102 is processed to provide a pluralityof features. For example including but not limited to HRV, complexanalysis, cardiac output data or the like. Most preferably featureextraction module 3 from the GPF vector for further analysis, forexample including the features depicted in Table 1 above

Machine learning module 4 is provided to classify the feature vector setprovided from module 3. Optionally, machine learning module 4 mayprovide for feature selection and dimensionality reduction wherein theGPF vector undergoes feature selection and dimensionality reduction toprovide a second vector therein providing a representation of the GPF ina smaller dimension. Optionally and preferably feature selection anddimensionality reduction techniques for example include but are notlimited to Multi Dimensional Scaling (MDS), Principal Component Analysis(PCA), Sparse PCA (SPCA), Fisher Linear Discriminant Analysis (FLDA),Sparse FLDA (SFLDA), Kernel PCA (KPCA), ISOMAP, Locally Linear Embedding(LLE), Laplacian Eigenmaps, Diffusion Maps, Hessian Eigenmaps,Independent Component Analysis (ICA), Factor analysis (FA), HierarchicalDimensionality Reduction (HDR), Sure Independence Screening (SIS),Fisher score ranks, t-test rank, Mann-Whitney U-test taken alone or inany combination thereof, or as known and accepted in the art.

Optionally classification is performed on at least one of the GPF vectoror a second vector comprising a dimensionally reduced version of GPFvector. Most preferably machine learning module 4 provides for painclassification into at least two or more classes, for example includingpain and non-pain groups. Optionally, and preferably a plurality ofoptional classifiers may be used for example including but not limitedNearest Shrunken Centroids (NSC), Classification and Regression Trees(CART), ID3, C4.5, Multivariate Additive regression splines (MARS),Multiple additive regression trees (MART), Nearest Centroid (NC),Shrunken Centroid Regularized Linear Discriminate and Analysis (SCRLDA),Random Forest, Boosting, Bagging Classifier, AdaBoost, RealAdaBoost,LPBoost, TotalBoost, BrownBoost, MadaBoost, LogitBoost, GentleBoost,RobustBoost, Support Vector Machine (SVM), kernelized SVM, Linearclassifier, Quadratic Discriminant Analysis (QDA) classifier, NaïveBayes Classifier and Generalized Likelihood Ratio Test (GLRT) classifierwith plug-in parametric or non-parametric class conditional densityestimation, k-nearest neighbor, Radial Base Function (RBF) classifier,Multilayer Perceptron classifier, Bayesian Network (BN) classifier,multi-class classifier adapted from binary classifier with one-vs-onemajority voting, one-vs-rest, Error Correcting Output Codes,hierarchical multi-class classification, Committee of classifiers or thelike as are known and accepted in the art or any combination thereof.

Optionally classification may be further provided with a graded scoringrelating the level of pain classified.

Output module 5 is provided to display or otherwise communicate theclassification results provided by machine learning module 4. Optionallyclassification results may be displayed in a plurality of formats forexample including printout, visual display cues, acoustic cues or thelike. Optionally results may be displayed in graded or class formats, orthe like.

Optionally output module may communicate the classification results toand external device for further processing, medical intervention or thelike, for example classification results may be communicated to a drugadministration device to automatically, semi-automatically, or manuallycontrol the drug delivery of a pain medicament, or to indicate to acaregiver to control, change, decrease or increase the dosage ordelivery of a pain reducing drug.

FIGS. 2A-C provide a further optional depiction of machine learningmodule 4 of the processing module 104 of FIG. 1 wherein the featurevector set provided from module 3 is processed, most preferablysequentially, with a plurality of machine learning sub-modules, forexample, preprocessing and normalization sub-module, feature selectionand dimensionality reduction sub-module and classification sub-module.Machine-learning sub-modules implements machine learning techniques andmethods as is known and accepted in the art. Optionally, sub-modules ofthe machine learning module comprising preprocessing and normalization,feature selection and dimensionality reduction may be reorganized andprovided in all possible combination.

FIG. 2B provides a depiction of a non limiting example of an optionalembodiment of FIG. 2A, also described in Example 1 below, wherein themachine learning module 4 of the processing module 104 of FIG. 1 isprovided by a collection of machine learning sub-modules. For example,preprocessing is provided by a plurality of actions including outlierremoval and signal normalization utilizing baseline zero mean unitvariance; feature selection is provided through the use of a Fisherscore; dimensionality reduction is provided by Sparse Fisher LinearDiscriminant Analysis (SFLDA) while classification is provided by aRealAdaboost classifier. Most preferably machine learning module 4 ofFIG. 1 comprises at least one preprocessing sub-module, at least onefeature selection and dimensionality reduction sub-module, and at leastone classification sub-module; of the optional techniques describedabove, in any combination thereof as is know and accepted in the art.

FIG. 2C provides a depiction of a non limiting example of an optionalembodiment of FIG. 2A, also described in Example 3 below, wherein themachine learning module 4 of the processing module 104 of FIG. 1 isprovided by a collection of machine learning sub-modules. For example,preprocessing is provided by a plurality of actions including outlierremoval and signal normalization utilizing baseline zero mean unitvariance; feature selection and dimensionality reduction is provided byHierarchical dimensionality reduction while classification is providedby a Random Forest classifier.

FIGS. 3A-E depicts optional hardware configurations of optional systemsaccording to optional embodiments of the present invention as describedin FIG. 1 above. FIGS. 3A-3E depicts different remote and localconfigurations of signal acquisition, processing and display. FIGS. 3Aand 3E depicts systems that may be adapted to provide for a fully localsystem, FIGS. 3B and 3C depicts systems that may be optionallyconfigured to be semi remote, while FIG. 3D depicts a system that may beadapted to be fully remote.

FIG. 3A depicts optional configuration system 300 comprising computer302 and an external signal acquisition and display monitor 308,optionally and preferably communicating using communication protocol310. Most preferably, monitor 308 comprises a stand alone or externalmonitor comprising signal acquisition and display that optionally andpreferably provides for signal acquisition that is preferably processedwith computer 302. Optionally and preferably computer 302 provides forthe processing for pain monitoring and detection according to thepresent invention. Optionally computer 302 may be realized as aprocessor, server or the like computing device as is known and acceptedin the art. Optionally communication protocol 310 mediates communicationand data exchange between computer 302 and external monitor 308providing for pain classification. Optionally, communication protocolmay be provided in a plurality of optional communication protocols as isknown and accepted in the art for example including but not limited towireless, wired, cellular, internet, Bluetooth, optical, IR or the likecommunication protocols as is known and accepted in the art. Optionallyand preferably computer 302 provides for pain classification andmonitoring according to the present invention while physiological signalacquisition and display is provided by monitor 308 most preferablymonitor 308 acquires the signals, transmit them to processor 302 andoptionally and preferably displays the classification result

FIG. 3C depicts an optional system 303 similar in configuration tosystem 300 of FIG. 3A above wherein the physiological signals arepreferably provided by both monitor 308 and acquisition transducers 304wherein both are communicated to computer 302 for processing accordingto an optional method according to the present invention.

FIG. 3B depicts optional configuration system 301 comprising computer302, signal acquisition transducers 304, physiological signalacquisition device 306, acquiring and displaying monitor 308 andcommunication protocol 310. Optionally signal acquisition transducer 304is most preferably provided in the form of optional sensors dependent onthe signal being acquired that is preferably coupled to signalacquisition device 306 most preferably to sample, amplify and processand record a physiological signal. Optionally and preferably device 306is a mobile, optionally provided in a plurality of formats for exampleincluding but not limited to mobile telephone, PDA, hand held device,MP3 player, dedicated or converted device adept for receiving andcommunicating a physiological signal. Optionally device 306 comprises aprocessor for performing initial signal processing of the physiologicalsignal. Optionally, device 306 does not comprise a processor adept forprocessing the sampled physiological signals and is linked usingcommunication protocol 310 to higher processing centers, for examplecomputer 302. Optionally, device 306 and computer 302 comprise a masterslave processing protocols for processing the physiological signalssampled with system acquisition transducers 304. Optionallycommunication protocol 310 may be further facilitated through aninternet connection. Optionally, computer 302 and display 308 may remoteto the signal acquisition transducer 304 and device 306 relaying oncommunication protocol 310 to communicate therebetween. Optionally,computer 302 and acquiring and/or displaying monitor 308 may beimplemented as a call center, telemedicine center, emergency medicalcenter, or remote center for remote pain classification and detectionfor pain management.

FIG. 3D depicts system 305 for remote pain classification comprising aremote signal acquisition device 306 is utilized to acquire andcommunicate physiological signals for pain classification provided forwith computer 302 provided in the form of a PDA comprising an intrinsicdisplay. System 305 therefore provides a remote system comprising remotesignal acquisition as well as remote signal processing and display.Conversely, FIG. 3E depicts system 307 comprising both local signalacquisition and processing.

FIG. 4 shows a flowchart of an exemplary method according to the presentinvention for pain monitoring and classification. In stage 1physiological signals are acquired, preferably a plurality of signalsare obtained as previously described in FIG. 1 using a plurality ofoptional sensors and/or transducers. Next in stage 2, signalpreprocessing is performed for example including but not limited to SNRoptimization, signal normalization, filtering or the like. Next in stage3, features extraction is performed to abstract and to preferably derivea great plurality of features GPF from the acquired physiologicalsignals. Next classification processed is initiated with a trainingprocess on the training set as depicted in stages 4′ and 5′. Mostpreferably an initial training process is performed in sage 4′ whereinfeature selection and dimensionality reduction functions are trained toimplement the optional dimensionality reduction and feature selectiontechniques as previously described. Next in stage 5′ the classifier istrained and implemented to identify the pain classification of thetraining set.

Optionally and preferably once the classifier has been trained and setfollowing stage 5′ future features obtained from the physiologicalsignals are classified without further training of the classifyingsystem of the present invention and implemented in stages 4 and 5. Instage 4 feature selection and dimensionality reduction is performedaccording to the dimensionality training determined in stage 4′. Next instage 5 classifications is performed with a classifier to determine thepain according to the dimensionality reduced GPF vector.

Finally in stage 6 pain monitoring and classification is displayedaccording to the appropriate time scale.

Example 1 Pain Classification Utilizing a Combination of FisherScore+SFLDA+RealAdaboost

The following is a non limiting example of the implementation of painmonitoring according to an optional embodiment of an optional method forpain monitoring and classification according to the present invention.The following example provides an illustrative example of painclassification for pain monitoring in a subject wherein the processingmethods comprise utilizing a priori data, feature extraction andselection based on Fisher Score rank, dimensionality reduction usingSFLDA, and classification using a RealAdaboost classifier within aBoosting framework.

Materials and Methods

An experiment was conducted in order to develop, validate performancesand evaluate efficacy of the classification method and system fornon-invasive automated pain monitoring according to the presentinvention. Primary outcome measures was to compare the pain monitoringresults with the subjective pain report measured by the visual NumericPain Scale (NPS) to a given pain stimulus with reports divided topain/no-pain (binary test).

The study included 26 healthy volunteers. Inclusion criteria were: (I)free from chronic pain of any type, (II) no medication usage except fororal contraceptives, (III) ability to understand the purpose andinstructions of the study, and (IV) blood pressure <140/90. Exclusioncriteria were: (I) any type of preexisting condition, (II) use ofmedications or recreational drugs, or (III) pregnancy. The study wasapproved by the local ethics committee, and a written informed consentwas obtained from all participants prior to the beginning of theexperiment.

The cold pressor test (CPT) and heat pain test were chosen to provokethe experimental pain. The cold pressor test apparatus (ChillSafe 8-30,ScanLaf A/S Denmark) is a temperature-controlled water bath with amaximum temperature variance of ±0.5° C., which is continuously stirredby a pump. Volunteers were asked to place their right foot (until abovethe ankle) in the CPT bath in a still position and maintain their footin the water for 1 min each session. A thermal testing analyzer (TSA)thermode of 30×30 mm (Medoc TSA-2001 device, Ramat Ishai, Israel) wasattached to the skin of the right forearm to initiate heat pain. Duringpain stimuli sessions the thermode was heated at 10° C./sec to thetarget temperature (39°-48.5° C.) and with a plateau lasting 60 sec.

In order to evaluate the volunteer's cold and heat pain sensitivity theywere exposed to a range of different temperatures and reported perceivedpain on a 0-100 numeric pain scale (NPS). For evaluation of subject'sheat pain sensitivity they were exposed to 11 heat stimuli, ranging from37° C. to 50° C. with increasing rate of 10° C./sec, each with a plateaulasting for 10 seconds. In order to evaluate subject's cold painsensitivity we expose subjects for one min to the CPT using watertemperature of 12° C. Subjects reported their pain every 10 sec. Thecold and heat apparatus were then appropriately calibrated to initiatefeeling of no-pain (NPS 0-30), and pain (NPS >70). Due to thelimitations of the GCP and ethics committee approval, theminimum/maximum temperatures used were 1° C. for the cold pain and 48.5°C. for the heat pain.

Subjects received a full explanation about the purpose and design of thestudy and signed a written informed consent form. Prior to the beginningof the experiment, the familiarizing and calibration sessions wereconducted. Next, the experimenter connected the sensors to the subjectand during the next 5 min physiological signals were rerecorded forbaseline normalization purposes. During the two sessions eachparticipant received 4 heat stimuli and 3 cold stimuli sessions, lastingfor 1 minute each with intervals of 10-15 minutes between stimuli and30-45 minutes between sessions. Volunteers were unaware of the stimuliintensity. The order of the heat pain stimuli was randomly assigned,while the order of the cold stimuli was progressive from low to highlevel of pain (to avoid adaptation). The “no-pain” stimuli (25° C. CPTand 39° C. heat sensor) was introduced with intention to stimulate asensory experience similar to pain sessions but without painfulstimulus.

The NPS report, the physiological signals and the intensity of stimuliwere recorded and synchronized for subsequent processing.

The physiological signals were recorded and stored in a personalcomputer by the BioPac MP 100 system (BioPac System Inc., CA, USA) andits companion software AcqKnowledge 3.9.1 (BioPac System Inc.). One-leadelectrocardiogram (ECG) signal, 2-channel electroencephalogram (EEG)signal from forehead, photoplethysmograph (PPG) signal from right handfinger, and one lead external electromyogram (EMG) signal from righttrapezius muscle were sampled with a frequency of 500 Hz. External skintemperature from the dorsum of the right hand, respiration, and galvanicskin response (GSR) from right hand fingers were sampled with afrequency of 32.5 Hz. In addition, continuous blood pressure signal wasnon-invasively measured using the Finometer MIDI (Finapres MedicalSystems BV, Amsterdam, The Netherlands) and data were recorded usingcompanion software BeatScope EASY (Finapres Medical Systems BV).Continuous blood pressure was sampled with a frequency of 200 Hz.

The recorded physiological signals were extracted, synchronized andprocessed in off-line way using Matlab®2009 scientific software (TheMathworks, Inc., MA, USA).

The data was processed offline on PC computer. All signals wereprocessed using routine signal processing methods for noise and artifactfiltering (Oppenhem & Shafer 1999). For some signals (EEG, ECG, etc.)were used signal-specific data processing methods (Rangayyan 2002,Sannei & Chambers 2007). All extracted parameters were averaged (ifapplicable) with non-overlapping windows of 10 sec. Features utilizedare summarized in Table 1 above were extracted from a plurality of rawphysiological data. The extracted features were deployed to training andvalidate classification algorithm.

The machine learning module implemented as described in greater detailin FIG. 2B above.

During the preprocessing and normalization stages the data was manuallyexamined by trained professional. Patients with over-sensitivity andunder-sensitivity were removed from training data. This exclusion fromtraining data set is a non-limiting example of a priori knowledge thatis incorporated into the system and method of the present invention mostpreferably to improve classification results. Nevertheless, the dataexcluded from training the optional classifier was used to test theclassifier.

Next, feature normalization was performed by removing the patient'sfeatures baseline mean and normalizing the patient's feature baselinevariability, in the following manner:

$\frac{X_{i} - {{avg}\left( X_{i}^{baseline} \right)}}{{std}\left( X_{i}^{baseline} \right)}$

Next Feature Selection step was performed. During the training phasefeatures with maximum Fisher score were selected. Fisher score of i′thfeature was defined as

$F_{i} = \frac{\left( {{{avg}\left( X_{i}^{pain} \right)} - {{avg}\left( X_{i}^{{no}\text{-}{pain}} \right)}} \right)^{2}}{{{var}\left( X_{i}^{pain} \right)} + {{var}\left( X_{i}^{{no}\text{-}{pain}} \right)}}$

The 100 best features with highest Fisher scores were predeterminedduring the training phase. During feature selection only the 100 highlyranked features were kept, the remaining features were removed fromconsideration.

Next, Dimensionality Reduction was performed. During the training phaseof the non-limiting and optional method according to the presentinvention a transformation matrix was determined based on the principleof Sparse Fisher Linear Discriminant Analysis (SFLDA). The methods forcalculation of a single SFLDA transformation vector were rigorouslydescribed in (Moghaddam, et al., 2006). wherein the transformationvector is calculated either by exhaustive search or using greedyalgorithms (forward or backward). Multiplication of the transformationvector by a vector of the selected features following feature selectionas described above, results in a number which represents the firstreduced discriminate dimension. For pain monitoring uniquespecification, novel iterative procedure for finding, multipletransformation vectors and discriminative components was used as follow:

-   -   Input: B—N×N pooled between-class covariance matrix, W—N×N        pooled within-class covariance matrix, k—required sparsity level        for transformation vectors, P—number of transformation vectors    -   Algorithm:    -   1. Set p=1, I^(p)=(1 . . . N)    -   2. Input B_(I) _(p) , W_(I) _(p) into SFLDA algorithm in order        to compute optimum sparsity pattern I_(opt) ^(p) sparse        transformation vector a_(opt) ^(p).    -   3. Set I^(p+1)=I^(p)\I_(opt) ^(p). Set p=p+1. If p<P go back to        step 2.    -   Output: Set of sparsity patterns {(I_(opt) ¹, . . . I_(opt)        ^(P)} and transformation matrix A=[a_(opt) ¹ . . . a_(opt)        ^(P)].        Transformation matrix A was multiplied by a vector of selected        features and transformed it into a new vector most preferably        having a significantly reduced dimension. This modification of        SFLDA is known as Sparse Linear Discriminant Component Analysis        (SLDCA). In an experiment each transformation vector had        sparsity level 5, and a total 10 such vectors were computed.        Therein according to this non-limiting example a feature set of        100 was reduced to 10 following an optional dimensionality        reducing step comprising SFLDA as described above

The final stage of classification according to an optional method of thepresent invention was set to classify the vector of reduced dimensionsinto a binary class of Pain and No-Pain classes. An optional classifieraccording to the present invention, RealAdaboost was chosen and trainedduring a training phase.

In order to assess performances of proposed method and apparatus theleave-one-out cross-validation (Hastie et al., 2009) scheme has beenused. In order to prevent the situation where the algorithm was bothtrained and validated on the same data (in our case same subject), thealgorithm was applied N times, where N is the number of subjects. Ateach run, data were included in the training set from all subjectsexcluding one, and then the trained algorithm was scored on data fromthis subject. The Test Error is estimated by averaging overclassification errors from each of the N runs.

The algorithm was also tested on patients, which were declared asoutliers. In such cases, the algorithm was trained on all non-outlierpatients and was tested on all outliers.

Results

The performance of the algorithm is presented in the Table 2. Theoverall agreement, sensitivity, and positive predictive values (PPV) arepresented, together with their respective 95% exact binomial confidenceintervals

TABLE 2 Results Pain level estimated pain no pain Total Pain level pain425 68 493 no pain 181 1528 1709 Total 606 1596

TABLE 3 Pain vs. No pain 95% exact binomial Percent confidence intervalOverall agreement 88.69% 87.29% 89.99% PPV pain 70.13% 66.31% 73.75% PPVno pain 95.74% 94.63% 96.68% Sensitivity pain 86.21% 82.84% 89.13%Sensitivity no pain 89.41% 87.85% 90.83%

Example 2 Feature Selection and Dimensionality Reduction

The following example is of a non-limiting example showing theimportance of machine learning techniques according to the presentinvention to improve and provide for pain monitoring, classification andidentification. During the clinical study described in Example 1 aboveall features described in Table 1 were measured. Physiological signalssuch as an accelerometer, SPO2 were excluded as the clinical study wasperformed on awake patients. Fisher scores were calculated to rank allavailable features, according to an optional method of the presentinvention.

Analysis of the features showing the following features as having a highFisher score, as summarized in Table 4 and FIGS. 5A-C.

As can be seen these parameters represent both an autonomic systemactivity of a patient (e.g. amplitudes of PPG signal, Count ofSpontaneous Fluctuations in GSR signal, H.R.V extracted from PPG or ECGand its spectral analysis) as well as behavioral activity, such asrespiration rate and deviation, muscle activity (which may suggest on acertain discomfort).

TABLE 4 Total Signal Temporal Feature Spectral Features Features GSRCount of Spontaneous Power of frequency 8 Fluctuations (Peaks), bandsVLF/LF/HF/ Weighed Sum of Peaks above HF (weighted by their amplitude),Signal Max value in window, Smoothed Signal Value in window PPG Mean andStandard Power of frequency 9 Deviation of amplitude bands HF and LF,and of Peak, Trough, and LF/HF ratio calculated Max Rate point. Meanfrom spectral analysis and Standard Deviation of Peak to Notch of AUC,Peak to Notch interval and Pulse interval, and Pulse Transition Time.Transition Time (Peak to R-ECG interval). Respiratory Mean and StandardPower of frequency Deviation of the bands VLF/LF/HF/ amplitude of Peak,Peak above HF to Peak interval (respiratory rate), Smoothed Signalvalue, Signal Max and Min values in Window Temperature Weighed Sum ofPeaks Power of frequency 9 (fluctuations) in bands VLF/LF/ window, Meanabove HF. amplitude of Peaks in window ECG Mean and Standard Power of HFand LF 4 Deviation of R peak bands from Heart rate variability EEG/EMGBarlow's mean Energy in band 30-70 9 amplitude, Hjorth's Hz, mean power,mobility, between Barlow's mean channels covariance (we frequency andspectral use 2-channel EEG), purity time differentials of a signal. EGGPower of below LF 1 band

For example, FIG. 5A depicted scatter plots of 5 random features showingtwo classes of pain namely pain (“o” open circles and) and non pain (“+”plus sign) for a classification attempted with the initial features thathave not undergone dimensionality reduction according to the presentinvention. Accordingly, inspection of the scatter plots clearly showsthat the two classes are not separable either by linearly ornon-linearly classification.

FIG. 5B depicts scatter plots of 5 features extracted from the GPFvector based on the highest Fisher scores and as depicted in Table 4:

-   -   1st Fisher feature is (combination): number of Peaks, Weighed        Sum of Peaks, Signal Max value in window, Smoothed Signal Value        in window of GSR signal.    -   2nd Fisher feature is: Power of HF band of GSR signal.    -   3rd Fisher feature is: Standard deviation of Peaks to peak        intervals of Respiratory signal.    -   4th Fisher feature is combination: Mean of Peak, Maximum Rate        Point, and Trough Heart Rate from PPG signal.    -   5th Fisher feature is: Standard deviation of amplitude of Peaks        of Respiratory signal.

Although inspection of FIG. 5B shows that scoring features providesimproved classification as can be seen by the reduced overlap betweenthe two classes depicted when compared to FIG. 5A however a substantialdegree of overlap is observed.

FIG. 5C depicts scatter plots of 5 SLDCA components whereindimensionality reduction technique utilized SLDCA substantially improvespain classification between the two classes.

1st SLDCA component comprises a combination of:

-   -   Standard deviation of Peaks to peak intervals of Respiratory        signal (3 Fisher feature)    -   Deviation of amplitude of Peaks of Respiratory signal (5 Fisher        feature)    -   Number of Peaks, Weighed Sum of Peaks, Signal Max value in        window, Smoothed Signal Value in window of GSR signal (1 Fisher        feature)    -   Power of LF band of PPG signal    -   Mean of Peak, Maximum Rate Point, and Trough, Heart Rate from        PPG signal (4 Fisher feature).

2nd SLDCA component comprises a combination of:

-   -   Mean of the amplitude of Peak of Respiratory signal    -   Standard Deviation of amplitude of Peak of Respiratory signal    -   Power of bellow LF band EGG    -   Between channels covariance EEG    -   ●Power of above HF band of GSR

3rd SLDCA component is a combination of:

-   -   First differential of EMG signal    -   Standard Deviation of Peak-to-Peak interval of EGG    -   Peak-to-Notch interval of PPG.    -   Power of HF band of peak-to-notch variability.    -   Power of VLF band Temperature

4th SLDCA component is a combination of:

-   -   Power of HF band of Respiratory signal.    -   Heart Rate from ECG signal    -   Standard Deviation of amplitude of Peak of EGG    -   Pulse transition time (Peak to R) from ECG and PPG    -   Standard Deviation of Peak-to-Notch interval of PPG

5th SLDCA component is a combination of:

-   -   Power of bellow VLF band of Respiratory    -   Power of LF band EGG    -   Standard Deviation of Pulse transition time (Peak to R) from ECG        and PPG    -   LF/HF ratio of PPG    -   Power of LF band Temperature

Although inspection of FIG. 5C shows that SLDCA dimensionality reductionprovides improved classification as can be seen by the reduced overlapbetween the two classes depicted when compared to FIG. 5B still asubstantial degree of overlap is observed, therefore classificationmodule is still required.

Example 3 Pain Classification Utilizing a Combination of HierarchicalDimensionality Reduction+Random Forest Classification

The following example is of a non-limiting example showing theimportance of machine learning techniques according to the presentinvention to improve and provide for pain monitoring, classification andidentification, as depicted in FIG. 2C wherein the system and method ofthe present invention comprises a wherein feature selection anddimensionality reduction is provided by Hierarchical DimensionalityReduction (HDR), and classification is provided by a Random Forestclassifier.

Subjects

The study included 36 healthy volunteers, 23 female and 13 male, aged 20to 38 years (mean (SD) 26(4.3)). Inclusion criteria were: (I) free fromchronic pain of any type, (II) no medication usage except for oralcontraceptives, (III) ability to understand the purpose and instructionsof the study, and (IV) blood pressure <140/90. Exclusion criteria were:(I) any type of preexisting condition, (II) use of medications orrecreational drugs, or (III) pregnancy. The study was approved by thelocal ethics committee, and a written informed consent was obtained fromall participants prior to the beginning of the experiment.

Instruments for Pain and Stress Stimulation

Assessment of Cold and Heat Pain Perception

The cold pressor test (CPT) and heat pain test were chosen to provokethe experimental pain. The cold pressor test apparatus (ChillSafe 8-30,ScanLaf A/S Denmark) is a temperature-controlled water bath with amaximum temperature variance of ±0.5° C., which is continuously stirredby a pump. Volunteers were asked to place their right foot (until abovethe ankle) in the CPT bath in a still position and maintain their footin the water for 1 min each session. A thermal testing analyzer (TSA)thermode of 30×30 mm (Medoc TSA-2001 device, Ramat Ishai, Israel) wasattached to the skin of the right forearm to initiate heat pain. Duringpain stimuli sessions the thermode was heated at 10° C./sec to thetarget temperature (39-48.5° C.) and with a plateau lasting 60 sec.

In order to evaluate the volunteer's cold and heat pain sensitivity theywere exposed to a range of different temperatures and reported perceivedpain on a 0-100 numeric pain scale (NPS). For evaluation of subject'sheat pain sensitivity they were exposed to 11 heat stimuli, ranging from37° C. to 50° C. with increasing rate of 10° C./sec, each with a plateaulasting for 10 seconds. In order to evaluate subject's cold painsensitivity we expose subjects for one min to the CPT using watertemperature of 12° C. Subjects reported their pain every 10 sec. Thecold and heat apparatus were then appropriately calibrated to initiatefeeling of no-pain (NPS 0-15), mild pain (NPS 15-45), moderate pain (NPS45-75), and severe pain (NPS >75). Due to the limitations of the GCP andethics committee approval, the minimum/maximum temperatures used were 1°C. for the cold pain and 48.5° C. for the heat pain. On severaloccurrences the calibrated temperatures of “moderate pain” weresignificantly close to minimum/maximum allowed temperatures. At suchcases “severe pain” stimuli temperatures were set to minimum/maximumtemperatures.

Mental Stress Protocol

The control sessions of mental stress was performed. In order to inducemental stress the mental arithmetic test Paced Auditory Serial AdditionTask (PASAT), (Gronwall & Wrightson 1974) has been used, which consistsof 1 min of auditory presentation of random digits from one to nine,with an interval of 2 sec. The subjects' were asked continuously expressthe sum of the two last digits.

Physiological Signals Recording

The physiological signals were recorded and stored in a personalcomputer by the BioPac MP 100 system (BioPac System Inc., CA, USA) andits companion software AcqKnowledge 3.9.1 (BioPac System Inc.). One-leadelectrocardiogram (ECG) signal, 2-channel electroencephalogram (EEG)signal from forehead, photoplethysmograph (PPG) signal from right handfinger, and one lead external electromyogram (EMG) signal from righttrapezius muscle were sampled with a frequency of 500 Hz. External skintemperature from the dorsum of the right hand, respiration, and galvanicskin response (GSR) from right hand fingers were sampled with afrequency of 32.5 Hz. In addition, continuous blood pressure signal wasnon-invasively measured using the Finometer MIDI (Finapres MedicalSystems BV, Amsterdam, The Netherlands) and data were recorded usingcompanion software BeatScope EASY (Finapres Medical Systems BV).Continuous blood pressure was sampled with a frequency of 200 Hz.

Data Processing and Parameters Extraction

The recorded physiological signals were extracted, synchronized andprocessed in off-line way using Matlab®2009 scientific software (TheMathworks, Inc., MA, USA).

All signals were processed using routine signal processing methods fornoise and artifact filtering (Oppenhem & Shafer 1999). For some signals(EEG, ECG, etc.) were used signal-specific data processing methods(Rangayyan 2002, Sannei & Chambers 2007). All extracted parameters wereaveraged (if applicable) with non-overlapping windows of 10 sec. Thedetailed list of parameters which were extracted from above mentionedphysiological signals can be found in Table 1

Normalization and Averaging

All continuous parameters were normalized by removing the parameter'sbaseline mean and normalizing the parameter's baseline variability, inthe following manner:

$\frac{X_{i} - {{avg}\left( X_{i}^{baseline} \right)}}{{std}\left( X_{i}^{baseline} \right)}$

In order to avoid bias and over-fitting, baseline signals of a specificsubject was not used for normalization of the subject's parameters. Inother words parameters of the subject where normalized based on baselineparameters of all other volunteers. Categorical parameters were notnormalized.

Finally, the parameters during each stimulus were averaged. Thus each 1min stimulus was represented by a vector of normalized parametersextracted from recorded physiological signals.

Dimensionality Reduction

All parameters with correlation r>0.8 were identified and consequentlyaveraged by Hierarchical Dimensionality Reduction (HDR) method Duda etal. 2000.

Classification Algorithm: Random Forest

The goal of the presented example was to demonstrate abilities of a painmonitoring device which optionally and preferably provides forpredicting the presence of pain sensation and classify its level. Afirst vector, according to optional embodiments of the presentinvention, comprising temporal, spectral, and other parameters extractedduring data processing and parameter extraction step constitutesprediction variables which we will use in prediction procedure.

In this example the classification algorithm which has been chosen wasthe Random Forest. Random Forest (Breiman 2001b) is a statisticallearning procedure that makes a prediction by aggregating the resultsfrom an ensemble of classification and regression trees (CART) (Breimanet al. 1984). Random Forest averages over multiple CART trees increasethe stability of the final algorithm. Each tree is grown on abootstrapped sample (random sampling with repetition) of originaltraining data. During the growing process of each tree, features arerandomly sampled as well to find a best split of each node of the tree.

Classification Algorithm Testing Methodology

In order to assess performances of proposed method and apparatus theleave-one-out cross-validation (Hastie et al., 2009) scheme has beenused. In order to prevent the situation where the algorithm was bothtrained and validated on the same data (in our case same subject), thealgorithm was applied N times, where N is the number of subjects. Ateach run, data were included in the training set from all subjectsexcluding one, and then the trained algorithm was scored on data fromthis subject. The Test Error is estimated by averaging overclassification errors from each of the N runs.

Study Design

Subjects received a full explanation about the purpose and design of thestudy and signed a written informed consent form. Prior to the beginningof the experiment, the familiarizing and calibration sessions wereconducted. Next, the experimenter connected the sensors to the subjectand during the next 5 min physiological signals were rerecorded forbaseline normalization purposes. After the baseline the physical stresssession was performed. During the two sessions each participant received3 heat stimuli and 3 cold stimuli sessions (mild pain, moderate pain andsevere pain), lasting for 1 minute each with intervals of 10-15 minutesbetween stimuli and 30-45 minutes between sessions. Volunteers wereunaware of the stimuli intensity. The order of the heat pain stimuli wasrandomly assigned, while the order of the cold stimuli was progressivefrom low to high level of pain (to avoid adaptation). The “no-pain”stimuli (25° C. CPT and 39° C. heat sensor) was introduced withintention to stimulate a sensory experience similar to pain sessions butwithout painful stimulus.

Results

The performance of the algorithm is assessed by presenting tables of thestimuli inflicted on the volunteers by experimenter versus the pain aspredicted and classified by the classification software based solely onrecorded physiological signals.

Since the cost function was used in the training of the classificationalgorithm, the weight matrices should be used to rank the estimationerror, i.e., to give a smaller penalty to a small error than a largeone; for example, in case of high pain stimulation, the deviceestimation of “no pain” will be considered a worse error than a “mildpain” estimation. If pain levels are ranked by their severity(increasing or decreasing order) the “weight” of estimated pain level igiven stimulated pain level j are defined according to the formula:w(i,j)=1−abs(i−j)/(N−1),where N is a total number of pain levels, and abs( ) denotes absolutevalue. This formula intends to introduce ranking constraint into painlevel estimation.

The weighting matrix used in the case of three categories is:

Pain level estimated Pain level stimulated Severe Mild No Pain Severe 1½ 0 Mild ½ 1 ½ No Pain 0 ½ 1

In case of two categories no weights are used, in other words all errorshad equal costs. In these cases the overall agreement, sensitivity, andpositive predictive values (PPV) are presented, together with theirrespective 95% exact binomial confidence intervals. Weights areintroduced in sensitivity and PPV calculations for multiclass cases.Thus, the resulted sensitivity and PPV values along with theircorresponded CI values are weighted sensitivity, PPV and CIrespectively. In these cases, the weighted overall agreement, weightedsensitivity and weighted PPV are presented, together with theirrespective 95% exact Wilson-score confidence intervals.

TABLE 5 Severe pain vs. No pain Pain level estimated severe pain no painTotal Pain level inflicted severe pain 249 33 282 no pain 28 271 299Total 277 304 581

TABLE 6 Severe pain vs. No pain 95% exact binomial Percent confidenceinterval Overall agreement 89.50% 86.72% 91.87% PPV severe pain 89.89%85.72% 93.18% PPV no pain 89.14% 85.09% 92.41% Sensitivity severe pain88.30% 83.96% 91.81% Sensitivity no pain 90.64% 86.75% 93.69%

TABLE 7 Severe pain vs. Mental Stress Pain level estimated severe painmental stress Total Pain level inflicted severe pain 292 37 329 mentalstress 69 260 329 Total 361 297 658

TABLE 8 Severe pain vs. Mental Stress 95% exact binomial Percentconfidence interval Overall agreement 83.89% 80.85% 86.62% PPV severepain 80.89% 76.44% 84.81% PPV mental stress 87.54% 83.24% 91.07%Sensitivity severe pain 88.75% 84.83% 91.96% Sensitivity mental stress79.03% 74.22% 83.30%

TABLE 9 Severe pain vs. Mild Pain vs. No pain Pain level estimatedsevere pain mild pain no pain Total Pain level severe pain 191 75 16 282inflicted mild pain 59 178 39 276 no pain 4 69 226 299 Total 254 322 281857

TABLE 10 Severe pain vs. Mild Pain vs. No pain 95% exact Wilson Percentconfidence interval Overall agreement 83.55% 80.92% 85.88% PPV severepain 86.81% 82.10% 90.43% PPV mild pain 77.64% 72.78% 81.85% PPV no pain87.37% 82.97% 90.75% Sensitivity severe pain 81.03% 76.05% 85.18%Sensitivity mild pain 82.25% 77.30% 86.30% Sensitivity no pain 87.12%82.85% 90.45%

These results show that the system and method implemented according tothe present invention is capable of classifying and detecting the painstatus in a patient as shown in Tables 6, 8 and 10

REFERENCES

-   Akselrod, S., D. Gordon, F. A. Ubel, D. C. Shannon, A. C. Berger,    and R. J. Cohen. “Power spectrum analysis of heart rate fluctuation:    a quantitative probe of beat-to-beat cardiovascular control.”    Science 213, no. 4504 (July 1981): 220-222.-   Belkin, M., and P. Niyogi. “Laplacian eigenmaps for dimensionality    reduction and data representation.” Neural Computation 15 (2003):    13731396.-   Bishop, C. M. Pattern recognition and machine learning. New York:    Springer, 2006.-   Borg, I., and P. Groenen. Modern Multidimensional Scaling: theory    and applications. New York: Springer-Verlag, 2005.-   Breiman L, Friedman J, Olshen R, Stone C. Classification and    Regression Trees. Belmont, Calif.: Wadsworth; 1984.-   Breiman, L. “Bagging Predictors.” Machine Learning 24, no. 2 (1996):    123-140.-   Breiman, L. “Random Forests.” Machine Learning V45, no. 1 (October    2001): 5-32.-   Breiman L. Random Forests, Machine Learning J., 2001b; 45:5-32.-   Bruhn, J., H. Ropcke, and A. Hoeft. “Approximate entropy as an    electroencephalographic measure of anesthetic drug effect during    desflurane anesthesia.” Anesthesiology 92, no. 3 (March 2000):    715-726.-   Chapelle, O., V. Vapnik, O. Bousquet, and S. Mukherjee. “Choosing    Kernel Parameters for Support Vector Machines.” Machine Learning,    2000.-   Cohen, D. H., and S. M. Sherman. “The autonomic nervous system and    its central control.” Edited by R. Berne, and M. Levy, 322. St.    Louis: C.V. Mosby, 1983.-   Coifman, R. R., et al. “Geometric diffusions as a tool for harmonic    analysis and structure definition of data.” Proceedings of the    National Academy of Sciences 102 (2005): 7426-7438.-   d'Aspremont, A., L. El Ghaoui, M. I. Jordan, and G. R. G. Lanckriet.    “A direct formulation for sparse PCA using semidefinite    programming.” Cambridge, Mass.: MIT Press, 2005.-   Deschamps, A., I. Kaufman, S. B. Backman, and G. Plourde. “Autonomic    nervous system response to epidural analgesia in laboring patients    by wavelet transform of heart rate and blood pressure variability.”    Anesthesiology 101, no. 1 (July 2004): 21-27.-   Duda R, Hart P, Stork D. Pattern Classification. New-York: Willey;    2000.-   Donoho, D. L., and C. Grimes. “Hessian eigenmaps: Locally linear    embedding techniques for high-dimensional data.” PNAS 100, no. 10    (May 2003): 5591-5596.-   Goncharova, I. I., and J. S. Barlow. “Changes in EEG mean frequency    and spectral purity during spontaneous alpha blocking.”    Electroencephalography and clinical neurophysiology 76, no. 3    (September 1990): 197-204.-   Gronwall D, Wrightson P. Delayed Recovery of intellectual function    after minor head injury. Lancet 2 (1974):605-609.-   Guignard, B. “Monitoring analgesia.” Best practice and research.    Clinical anaesthesiology 20, no. 1 (March 2006): 161-180.-   Guo, Y., T. Hastie, and R. Tibshirani. “Regularized linear    discriminant analysis and its application in microarrays.”    Biostatistics (Oxford University Press) 8, no. 1 (January 2007):    86-100.-   Guyton, A. C. Human Physiology and Mechanisms of Disease.    Philadelphia: Saunders, 1982.-   Fagrell, B. Microcirculation of the skin. The physiology and    pharmacology of the microcirculation. Vol. 2. New York: Academic    Press, 1984.-   Hastie, T., R. Tibshirani, and J. Friedman. The elements of    statistical learning: data mining, inference, and prediction. New    York: Springer, 2009.-   Hjorth, B. “The physical significance of time domain descriptors in    EEG analysis.” Electroencephalography and clinical neurophysiology    34, no. 3 (March 1973): 321-325.-   Kohavi, R., and G. H. John. “Wrappers for Feature Subset Selection.”    Artificial Intelligence 97, no. 1-2 (1997): 273324.-   Lidberg, L., and G. B. Wallin. “Sympathetic Skin Nerve Discharges in    Relation to Amplitude of Skin Resistance Responses.”    Psychophysiology 18, no. 3 (1981): 268-270.-   Malik, M. “Heart rate variability: standards of measurement,    physiological interpretation and clinical use. Task Force of the    European Society of Cardiology and the North American Society of    Pacing and Electrophysiology.” Circulation 93, no. 5 (March 1996):    1043-1065.-   Moghaddam, B., Y. Weiss, and S. Avidan. “Spectral Bounds for Sparse    PCA: Exact and Greedy Algorithms.” Edited by Weiss, Y.,    Scholkopf, B. and Platt, J., 915-922. Cambridge, Mass.: MIT Press,    2006.-   Moghaddam, B., Y. Weiss, and S. Avidan. “Generalized Spectral Bounds    for Sparse LDA.” ICML. 2006a. 641-648.-   Oppenheim A V, Schafer R W. Discrete-time signal Processing. Upper    Saddle River, N.J.: Prentice Hall; 1999.-   Pagani, M., O. Rimoldi, and A. Malliani. “Low-frequency components    of cardiovascular variabilities as markers of sympathetic    modulation.” Trends in pharmacological sciences 13, no. 2 (February    1992): 50-54.-   Pan, J., and W. J. Tompkins. “A real-time QRS detection algorithm.”    IEEE Trans Biomed Eng 32, no. 3 (March 1985): 230-236.-   Pop-Jordanova, N., and J. Pop-Jordanov. “Spectrum-weighted EEG    frequency (“brain-rate”) as a quantitative indicator of mental    arousal.” Prilozi Makedonska akademija na naukite i umetnostite,    Oddelenie za biolovski i medicinski nauki 26, no. 2 (December 2005):    35-42.-   Rangayyan R M. Biomedical signal analysis: A case study approach.    IEEE Press, 2002.-   Sanei S, Chambers J A. EEG signal processing. New York: John Wiley &    Sons; 2007.-   Schlogl, Alois. “A comparison of multivariate autoregressive    estimators.” Signal Process. (Elsevier North-Holland, Inc.) 86, no.    9 (September 2006): 2426-2429-   Scholkopf, B., A. Smola, and K.-R. Muller. “Nonlinear component    analysis as a kernel eigenvalue problem.” Neural Comput. (MIT Press)    10, no. 5 (1998): 1299-1319.-   Seitsonen, E. R., I. K. Korhonen, M. J. Van Gils, J. M.    Lotjonen, K. T. Kortilla, and A. M. Yli-Hankala. “EEG spectral    entropy, heart rate, photoplethysmography and motor responses to    skin incision during sevoflurane anaesthesia.” Acta Anaesthesiol    Scand, 49, (2005): 284292.-   Roweis, S., and L. Saul. “Nonlinear dimensionality reduction by    locally linear embedding.” Science 290. (2000): 2323-2326.-   Tenenbaum, J., V. de Silva, and J. C. Langford. “A global geometric    framework for nonlinear dimensionality reduction.” Science 290    (2000): 23192323.-   Thomas, D. W., and J. M. Evans. “Lower oesophageal contractility    monitoring during anaesthesia for cardiac surgery: preliminary    observations.” Annals of the Royal College of Surgeons of England    71, no. 5 (September 1989): 311-315.-   Tibshirani, R., T. Hastie, B. Narasimhan, and G. Chu. “Diagnosis of    multiple cancer types by shrunken centroids of gene expression.”    PNAS 99, no. 10 (May 2002): 6567-6572.-   van den Berg, A. A., D. Savva, and N. M. Honjol. “Attenuation of the    haemodynamic responses to noxious stimuli in patients undergoing    cataract surgery. A comparison of magnesium sulphate, esmolol,    lignocaine, nitroglycerine and placebo given i.v. with induction of    anaesthesia.” European Journal of Anaesthesiology 14, no. 02 (2006):    134-147.-   Vapnik, V. N. Statistical Learning Theory. Edited by Haykin S. New    York: Wiley, 1998.-   Wackermann, J. “Towards a quantitative characterisation of    functional states of the brain: from the non-linear methodology to    the global linear description.” International Journal of    Psychophysiology 34, no. 1 (October 1999): 65-80.-   Weiss, T., A. Del Bo, N. Reichek, and K. Engelman. “Pulse transit    time in the analysis of autonomic nervous system effects on the    cardiovascular system.” Psychophysiology 17, no. 2 (March 1980):    202-207.-   Zou, H., T. Hastie, and R. Tibshirani. “Sparse Principal Component    Analysis.” Journal of Computational Graphical Statistics 15 (2006):    265286.

While the invention has been described with respect to a limited numberof embodiments, it will be appreciated that many variations,modifications and other applications of the invention may be made.

What is claimed is:
 1. A method for monitoring pain of a patient, themethod comprising: obtaining at least two physiological signalscomprising blood volume change and at least one signal selected from thegroup consisting of: Galvanic Skin Response (GSR); electrocardiogram(ECG), respiration, internal body temperature, skin temperature,electrooculography (EOG), pupil diameter, electroencephalogram (EEG),frontalis electromyogram (FEMG), electromyography (EMG),electro-gastro-gram (EGG), partial pressure of carbon dioxide, andaccelerometer readings; processing said at least two physiologicalsignals to improve signal quality, thereby obtaining at least twoprocessed signals; generating a first vector, said first vectorcomprising at least three features extracted from said at least twophysiological signals, wherein said at least three features comprise atleast one feature selected from the group consisting of: mean Peak (P)amplitude, Peak (P) mean amplitude, Peak (P) std of amplitude, Through(T) amplitude, Trough (T) mean amplitude, Trough (T) std of amplitude,peak to peak intervals, Peak-to-Peak High Frequency (P-P HF) Power, andany combination thereof; transforming said first vector into a secondvector, said transformation comprising normalization; and monitoring thepain status of the patient by applying a classification algorithmadapted to classify said second vector into a graduated scalerepresenting the level of pain; wherein said classification algorithmcomprises an ensemble of classification and regression trees.
 2. Themethod of claim 1, wherein said pain monitoring further comprisescommunicating said monitored pain to a receiving unit selected from thegroup consisting of: a higher processing center, person, caregiver, callcenter and any combination thereof.
 3. The method of claim 1, furthercomprising obtaining and processing a priori data, wherein said a prioridata is selected from the group consisting of: environmental parameters,patient parameters, disease, stimulus, medicament and any combinationthereof.
 4. The method according to claim 1, wherein said ensemble ofclassification and regression trees, comprises a random forestclassifier or a boosting framework.
 5. The method of claim 1, whereinsaid classification algorithm is adapted for pain experienced with aparticular disease, stimulus or medicament.
 6. The method of claim 1,wherein said patient is in a state of consciousness selected from thegroup consisting of: unconscious, under general anesthesia, sedated,partially sedated, awake, and semi-awake.
 7. A system for monitoring apain of a patient, said system comprising: a signal acquisition modulecomprising at least one sensor and/or transducers for measuring and/orobtaining at least two physiological signals comprising blood volumechange and at least one signal selected from the group consisting of:Galvanic Skin Response (GSR); electrocardiogram (ECG), respiration,internal body temperature, skin temperature, electrooculography (EOG),pupil diameter, electroencephalogram (EEG), frontalis electromyogram(FEMG), electromyography (EMG), electro-gastro-gram (EGG), partialpressure of carbon dioxide, and accelerometer readings; and a processingmodule for processing said at least two physiological signals, saidprocessing comprising: i. processing said at least two physiologicalsignal to improve signal quality, thereby forming at least two processedsignals; ii. generating a first vector, said first vector comprising atleast three features extracted from said at least two physiologicalsignals, wherein said at least three features comprise at least onefeature selected from the group consisting of: mean Peak (P) amplitude,Peak (P) mean amplitude, Peak (P) std of amplitude, Through (T)amplitude, Trough (T) mean amplitude, Trough (T) std of amplitude, peakto peak intervals, Peak-to-Peak High Frequency (P-P HF) Power, and anycombination thereof; iii. transforming said first vector into a secondvector said transformation comprising normalization; and iv. monitoringthe pain of the patient by applying a classification algorithm adaptedto classify said second vector into graduated scale representing thelevel of pain; wherein said classification algorithm comprises anensemble of classification and regression trees.
 8. The system of claim7, further comprising a display module for displaying said monitoredpain.
 9. The method of claim 1, wherein said at least three featuresfurther comprise PPG maximum rate (M) point, PPG dicrotic notch (N), PPGPP/PT/PN/NT/NM intervals, mean and std (variability) of intervals, PPGPP variability (PPG-HRV) in frequency bands: VLF, LF, MF, HF and LF/HF,area under the PPG curve (AUC), PPG spectrum envelope, PPG-HRV waveletanalysis, PPG-RSA (respiratory sinus arrhythmia), GSR amplitude, meanamplitude and std (variability) of amplitude, GSR PP interval, mean andstd (variability) of interval, GSR phasic EDA amplitude, mean amplitudeand std (variability) of amplitude, spectrum of GSR signal, peakfrequency, GSR wavelet analysis, ECG Q/R/S/T/P amplitude mean and std(variability) of amplitude, ECG RR/PQ/PR/QT/RS/ST interval, mean and std(variability) of interval, ECG RR variability (ECG-HRV) in frequencybands: VLF, LF, MF, HF and LF/HF, ECG-RSA (respiratory sinusarrhythmia), ECG-HRV wavelet analysis, ECG-PPG pulse transition time(PTT), temperature amplitude, mean and std of amplitude, temperature PPinterval, mean and std (variability) of interval, temperature spectrum,temperature peak frequency, upper temperature peak amplitude, meanamplitude and std (variability) of amplitude, lower temperature peakamplitude mean amplitude and std (variability) of amplitude, respiratoryrate, spectrum analysis of respiratory signal, mean rate and std(variability) of rate, EEG/EMG power of: alpha, beta, gamma, delta andtheta frequency bands, EEG/EMG mean frequency, EEG/EMG peak frequency,EEG/EMG total power, EEG/EMG spectral edge frequency, EEG/EMGapproximate entropy, EEG/EMG BSR (burst suppression ratio), EEG/EMGBcSEF, EEG/EMG WSMF, EEG/EMG CUP, EEG/EMG SpEn, EEG/EMG BcSpEn, EEG/EMGbeta ratio, EEG/EMG histogram parameters, EEG/EMG AR parameters, EEG/EMGnormalized slope description (Hjorth parameters), EEG/EMG Barlowparameters, EEG/EMG Wackerman parameters, EEG/EMG brain rate, EEG/EMGSynchFastSlow, EEG/EMG OMT EEG/EMG spectrum analysis, BP Spectrumanalysis EMG SLOC, average/variability of end-tidal airway gasses,average of accelerometer X, Y, Z, theta, accelerometer movementanalysis, coherence between 2 or more EEG/FEMG channels and combinationsthereof.
 10. The method of claim 1, wherein said at least three featuresfurther comprise GSR amplitude, GSR PP time intervals, temperatureamplitude, temperature peak (P) amplitude, temperature PP timeintervals, ECG-PPG PTT, ECG RR time intervals, ECG-HRV power of VLF, LFand HF frequency bands, respiration rate, spectrum analysis ofrespiratory signal, EEG/EMG spectrum analysis, BP spectrum analysis,coherence between 2 or more EEG/FEMG channels, accelerometer movementanalysis; and combinations thereof.
 11. The system of claim 7, whereinsaid at least three features further comprise PPG maximum rate (M)point, PPG dicrotic notch (N), PPG PP/PT/PN/NT/NM intervals, mean andstd (variability) of intervals, PPG PP variability (PPG-HRV) infrequency bands: VLF, LF, MF, HF and LF/HF, area under the PPG curve(AUC), PPG spectrum envelope, PPG-HRV wavelet analysis, PPG-RSA(respiratory sinus arrhythmia), GSR amplitude, mean amplitude and std(variability) of amplitude, GSR PP interval, mean and std (variability)of interval, GSR phasic EDA amplitude, mean amplitude and std(variability) of amplitude, spectrum of GSR signal, peak frequency, GSRwavelet analysis, ECG Q/R/S/T/P amplitude mean and std (variability) ofamplitude, ECG RR/PQ/PR/QT/RS/ST interval, mean and std (variability) ofinterval, ECG RR variability (ECG-HRV) in frequency bands: VLF, LF, MF,HF and LF/HF, ECG-RSA (respiratory sinus arrhythmia), ECG-HRV waveletanalysis, ECG-PPG pulse transition time (PTT), temperature amplitude,mean and std of amplitude, temperature PP interval, mean and std(variability) of interval, temperature spectrum, temperature peakfrequency, upper temperature peak amplitude, mean amplitude and std(variability) of amplitude, lower temperature peak amplitude meanamplitude and std (variability) of amplitude, respiratory rate, spectrumanalysis of respiratory signal, mean rate and std (variability) of rate,EEG/EMG power of: alpha, beta, gamma, delta and theta frequency bands,EEG/EMG mean frequency, EEG/EMG peak frequency, EEG/EMG total power,EEG/EMG spectral edge frequency, EEG/EMG approximate entropy, EEG/EMGBSR (burst suppression ratio), EEG/EMG BcSEF, EEG/EMG WSMF, EEG/EMG CUP,EEG/EMG SpEn, EEG/EMG BcSpEn, EEG/EMG beta ratio, EEG/EMG histogramparameters, EEG/EMG AR parameters, EEG/EMG normalized slope description(Hjorth parameters), EEG/EMG Barlow parameters, EEG/EMG Wackermanparameters, EEG/EMG brain rate, EEG/EMG SynchFastSlow, EEG/EMG OMTEEG/EMG spectrum analysis, BP Spectrum analysis EMG SLOC,average/variability of end-tidal airway gasses, average of accelerometerX, Y, Z, theta, accelerometer movement analysis, coherence between 2 ormore EEG/FEMG channels and combinations thereof.
 12. The system of claim7, wherein said at least three features further comprise GSR amplitude,GSR PP time intervals, temperature amplitude, temperature peak (P)amplitude, temperature PP time intervals, ECG-PPG PTT, ECG RR timeintervals, ECG-HRV power of VLF, LF and HF frequency bands, respirationrate, spectrum analysis of respiratory signal, EEG/EMG spectrumanalysis, BP spectrum analysis, coherence between 2 or more EEG/FEMGchannels, accelerometer movement analysis; and combinations thereof. 13.The system according to claim 7, wherein said ensemble of classificationand regression trees comprises a random forest classifier or a boostingframework.
 14. The system of claim 7, wherein said signal acquisitionmodule is further adapted to obtain a priori data selected from thegroup consisting of: environmental parameters, patient parameters,disease, stimulus, medicament and any combination thereof, and whereinsaid processing module is further adapted to process said priori data.15. The system of claim 7, further comprising a communication module forcommunicating said monitored pain of the patient to a receiving unitselected from the group consisting of a higher processing center, aperson, a caregiver, a call center and any combination thereof.